Соберите матрицы конечного элемента
FEM = assembleFEMatrices(___,state)state массив структур. Функция использует time поле структуры для зависящих от времени моделей и поле u решения для нелинейных моделей. Можно использовать этот аргумент с любым из предыдущих синтаксисов.
Создайте модель PDE для уравнения Poisson на L-образной мембране с нулем граничные условия Дирихле.
model = createpde(1); geometryFromEdges(model,@lshapeg); specifyCoefficients(model,'m',0,'d',0,'c',1,'a',0,'f',1); applyBoundaryCondition(model,'Edge',1:model.Geometry.NumEdges, ... 'u',0);
Сгенерируйте mesh и получите матрицы конечного элемента по умолчанию для проблемы и mesh.
generateMesh(model,'Hmax',0.2);
FEM = assembleFEMatrices(model)FEM = struct with fields:
K: [401x401 double]
A: [401x401 double]
F: [401x1 double]
Q: [401x401 double]
G: [401x1 double]
H: [80x401 double]
R: [80x1 double]
M: [401x401 double]
Сделайте расчеты быстрее путем определения который матрицы конечного элемента собраться.
Создайте переходную тепловую модель и включайте геометрию встроенной функции squareg.
thermalmodel = createpde('thermal','steadystate'); geometryFromEdges(thermalmodel,@squareg);
Постройте геометрию с метками ребра.
pdegplot(thermalmodel,'EdgeLabels','on') xlim([-1.1 1.1]) ylim([-1.1 1.1])
Задайте теплопроводность материала и внутреннего источника тепла.
thermalProperties(thermalmodel,'ThermalConductivity',0.2);
internalHeatSource(thermalmodel,10);Установите граничные условия.
thermalBC(thermalmodel,'Edge',[1,3],'Temperature',100);
Сгенерируйте mesh.
generateMesh(thermalmodel);
Соберите жесткость и большие матрицы.
FEM_KM = assembleFEMatrices(thermalmodel,'KM')FEM_KM = struct with fields:
K: [1541x1541 double]
M: [1541x1541 double]
Теперь соберите матрицы конечного элемента M, K, A, и F.
FEM_MKAF = assembleFEMatrices(thermalmodel,'MKAF')FEM_MKAF = struct with fields:
M: [1541x1541 double]
K: [1541x1541 double]
A: [1541x1541 double]
F: [1541x1 double]
Эти четыре матрицы M, K, A, и F соответствуют дискретизированным версиям коэффициентов УЧП m, c, a, и f. Эти четыре матрицы также представляют область модели конечного элемента УЧП. Вместо того, чтобы задать их явным образом, можно использовать domain аргумент.
FEMd = assembleFEMatrices(thermalmodel,'domain')FEMd = struct with fields:
M: [1541x1541 double]
K: [1541x1541 double]
A: [1541x1541 double]
F: [1541x1 double]
Эти четыре матрицы Q, G, H, и R, соответствуют дискретизированным версиям q, g, h, и r в спецификации граничного условия Неймана и Дирихле. Эти четыре матрицы также представляют контур модели конечного элемента УЧП. Используйте boundary аргумент, чтобы собрать только эти матрицы.
FEMb = assembleFEMatrices(thermalmodel,'boundary')FEMb = struct with fields:
H: [74x1541 double]
R: [74x1 double]
G: [1541x1 double]
Q: [1541x1541 double]
nullspace и stiff-spring МетодыСоздайте модель PDE для уравнения Poisson на L-образной мембране с нулем граничные условия Дирихле.
model = createpde(1); geometryFromEdges(model,@lshapeg); specifyCoefficients(model,'m',0,'d',0,'c',1,'a',0,'f',1); applyBoundaryCondition(model,'Edge',1:model.Geometry.NumEdges, ... 'u',0);
Сгенерируйте mesh.
generateMesh(model,'Hmax',0.2);Получите матрицы конечного элемента после наложения граничного условия с помощью подхода пустого пробела. Этот подход устраняет степени свободы Дирихле и обеспечивает уменьшаемую систему уравнений.
FEMn = assembleFEMatrices(model,'nullspace')FEMn = struct with fields:
Kc: [321x321 double]
Fc: [321x1 double]
B: [401x321 double]
ud: [401x1 double]
M: [321x321 double]
Получите решение УЧП с помощью nullspace матрицы конечного элемента.
un = FEMn.B*(FEMn.Kc\FEMn.Fc) + FEMn.ud;
Сравните этот результат с решением, данным solvepde. Эти два решения идентичны.
u1 = solvepde(model); norm(un - u1.NodalSolution)
ans = 0
Получите матрицы конечного элемента после наложения граничного условия с помощью жестко-пружинного подхода. Этот подход сохраняет степени свободы Дирихле, но налагает большой штраф на них.
FEMs = assembleFEMatrices(model,'stiff-spring')FEMs = struct with fields:
Ks: [401x401 double]
Fs: [401x1 double]
M: [401x401 double]
Получите решение УЧП с помощью жестко-пружинных матриц конечного элемента. Этот метод дает менее точное решение.
us = FEMs.Ks\FEMs.Fs; norm(us - u1.NodalSolution)
ans = 0.0098
Соберите матрицы конечного элемента для первых и последних временных шагов переходной структурной проблемы.
Создайте переходную структурную модель для решения твердой (3-D) проблемы.
structuralmodel = createpde('structural','transient-solid');
Создайте геометрию и включайте ее в модель. Постройте геометрию.
gm = multicylinder(0.01,0.05); addVertex(gm,'Coordinates',[0,0,0.05]); structuralmodel.Geometry = gm; pdegplot(structuralmodel,'FaceLabels','on','FaceAlpha',0.5)
Задайте модуль Молодежи и отношение Пуассона.
structuralProperties(structuralmodel,'Cell',1,'YoungsModulus',201E9, ... 'PoissonsRatio',0.3, ... 'MassDensity',7800);
Укажите, что нижняя часть цилиндра является фиксированным контуром.
structuralBC(structuralmodel,'Face',1,'Constraint','fixed');
Задайте гармоническое давление на верхнюю часть цилиндра.
structuralBoundaryLoad(structuralmodel,'Face',2,... 'Pressure',5E7, ... 'Frequency',50);
Определите нулевое начальное перемещение и скорость.
structuralIC(structuralmodel,'Displacement',[0;0;0], ... 'Velocity',[0;0;0]);
Сгенерируйте линейную mesh.
generateMesh(structuralmodel,'GeometricOrder','linear'); tlist = linspace(0,1,300);
Соберите матрицы конечного элемента для начального временного шага.
state.time = tlist(1); FEM_domain = assembleFEMatrices(structuralmodel,state)
FEM_domain = struct with fields:
K: [6609x6609 double]
A: [6609x6609 double]
F: [6609x1 double]
Q: [6609x6609 double]
G: [6609x1 double]
H: [252x6609 double]
R: [252x1 double]
M: [6609x6609 double]
Давление, поданное наверху цилиндра, является единственным зависящим от времени количеством в модели. Чтобы смоделировать динамику системы, соберите матрицу G конечного элемента граничной загрузки для начальной буквы, промежуточного звена и итоговых временных шагов.
state.time = tlist(1);
FEM_boundary_init = assembleFEMatrices(structuralmodel,'G',state)FEM_boundary_init = struct with fields:
G: [6609x1 double]
state.time = tlist(floor(length(tlist)/2));
FEM_boundary_half = assembleFEMatrices(structuralmodel,'G',state)FEM_boundary_half = struct with fields:
G: [6609x1 double]
state.time = tlist(end);
FEM_boundary_final = assembleFEMatrices(structuralmodel,'G',state)FEM_boundary_final = struct with fields:
G: [6609x1 double]
Соберите матрицы конечного элемента для проблемы теплопередачи с температурно-зависимой теплопроводностью.
Создайте установившуюся тепловую модель.
thermalmodelS = createpde('thermal','steadystate');
Создайте 2D геометрию путем рисования одного прямоугольника размер блока и второго прямоугольника размер паза.
r1 = [3 4 -.5 .5 .5 -.5 -.8 -.8 .8 .8]; r2 = [3 4 -.05 .05 .05 -.05 -.4 -.4 .4 .4]; gdm = [r1; r2]';
Вычтите второй прямоугольник сначала, чтобы создать блок с пазом.
g = decsg(gdm,'R1-R2',['R1'; 'R2']');
Преобразуйте decsg формат в геометрический объект. Включайте геометрию в модель и постройте геометрию.
geometryFromEdges(thermalmodelS,g); figure pdegplot(thermalmodelS,'EdgeLabels','on'); axis([-.9 .9 -.9 .9]);
Установите температуру на левом крае до 100 градусов. Установите поток тепла из блока на правом краю к-10. Верхние и нижние ребра и ребра в полости все изолируются: нет никакой теплопередачи через эти ребра.
thermalBC(thermalmodelS,'Edge',6,'Temperature',100); thermalBC(thermalmodelS,'Edge',1,'HeatFlux',-10);
Задайте теплопроводность материала как простая линейная функция температурного u.
k = @(~,state) 0.7+0.003*state.u;
thermalProperties(thermalmodelS,'ThermalConductivity',k);Сгенерируйте mesh.
generateMesh(thermalmodelS);
Вычислите установившееся решение.
Rnonlin = solve(thermalmodelS);
Поскольку теплопроводность нелинейна (зависит от температуры), вычислите системные матрицы, соответствующие сходившейся температуре. Присвойте температурное распределение u поле state массив структур. Поскольку u поле должно содержать вектор-строку, транспонировать температурное распределение.
state.u = Rnonlin.Temperature.';
Соберите матрицы конечного элемента с помощью температурного распределения в узлах.
FEM = assembleFEMatrices(thermalmodelS,'nullspace',state)FEM = struct with fields:
Kc: [1277x1277 double]
Fc: [1277x1 double]
B: [1320x1277 double]
ud: [1320x1 double]
M: [1277x1277 double]
Вычислите решение с помощью системных матриц, чтобы проверить, что они дают к той же температуре как Rnonlin.
u = FEM.B*(FEM.Kc\FEM.Fc) + FEM.ud;
Сравните этот результат с решением, данным solve.
norm(u - Rnonlin.Temperature)
ans = 7.2809e-05
PDEModel | ThermalModel | StructuralModel | ElectromagneticModel | solvepde | solve