Собрать матрицы конечных элементов
FEM = assembleFEMatrices(___,state)state структурный массив. Функция использует time поле структуры для зависящих от времени моделей и поле решения u для нелинейных моделей. Этот аргумент можно использовать с любым из предыдущих синтаксисов.
Создайте модель PDE для уравнения Пуассона на 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);
Создайте сетку и получите матрицы конечных элементов по умолчанию для задачи и сетки.
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]
T: [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);
Создайте сетку.
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 соответствуют дискретизированным версиям коэффициентов PDE m, c, a и F. Эти четыре матрицы также представляют область конечноэлементной модели PDE. Вместо явного указания их можно использовать 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 в спецификации граничных условий Неймана и Дирихле. Эти четыре матрицы также представляют границу конечноэлементной модели PDE. Используйте 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 для уравнения Пуассона на 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);
Создайте сетку.
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]
Получение решения для PDE с помощью 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]
Получите решение для PDE, используя матрицы конечных элементов жесткой пружины. Этот метод дает менее точное решение.
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]);
Создайте линейную сетку.
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');
Создайте 2-D геометрию, нарисовав один прямоугольник размером блока и второй прямоугольник размером паза.
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);Создайте сетку.
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 = 9.8994e-05
ElectromagneticModel | PDEModel | solve | solvepde | StructuralModel | ThermalModel