Установите начальные условия для модели с мелкой сеткой при помощи решения крупной сетки от предыдущего анализа.
Создайте модель PDE и включайте геометрию встроенной функции squareg
.
model = createpde; geometryFromEdges(model,@squareg);
Задайте коэффициенты, примените граничные условия и установите начальные условия.
specifyCoefficients(model,'m',0,'d',1,'c',5,'a',0,'f',0.1); applyBoundaryCondition(model,'dirichlet','Edge',1,'u',1); setInitialConditions(model,10);
Сгенерируйте сравнительно крупную сетку с целевой максимальной длиной ребра элемента 0,1.
generateMesh(model,'Hmax',0.1);
Решите модель для целого отрезка времени 0 в течение 0,02 секунд.
tlist = linspace(0,2E-2,20); Rtotal = solvepde(model,tlist);
Интерполируйте решение в начале координат для целого отрезка времени.
singleSpanSol = Rtotal.interpolateSolution(0,0,1:numel(tlist));
Теперь решите модель для первой половины отрезка времени. Вы будете использовать это решение в качестве начального условия при решении модели с более прекрасной mesh для второй половины отрезка времени.
tlist1 = linspace(0,1E-2,10); R1 = solvepde(model,tlist1);
Создайте interpolant, чтобы интерполировать начальное условие.
x = model.Mesh.Nodes(1,:)'; y = model.Mesh.Nodes(2,:)'; interpolant = scatteredInterpolant(x,y,R1.NodalSolution(:,end));
Сгенерируйте более прекрасную mesh путем установки целевой максимальной длины ребра элемента на 0,05.
generateMesh(model,'Hmax',0.05);
Используйте результаты модели крупной сетки в качестве начального условия для модели с более прекрасной mesh. Для определения icFcn
функционируйте, смотрите Функцию Начальных условий.
setInitialConditions(model,@(region) icFcn(region,interpolant));
Решите модель для второй половины отрезка времени.
tlist2 = linspace(1E-2,2E-2,10); R2 = solvepde(model,tlist2);
Интерполируйте решения в начале координат для первого и вторых половин отрезка времени.
multispanSol1 = R1.interpolateSolution(0,0,1:numel(tlist1)); multispanSol2 = R2.interpolateSolution(0,0,1:numel(tlist2));
Постройте все три решения в начале координат.
figure plot(tlist,singleSpanSol) hold on plot(tlist1, multispanSol1,'r*') plot(tlist2, multispanSol2,'ko') legend('Overall solution','Coarse mesh solution', 'Fine mesh solution')
function u0 = icFcn(region,interpolant) u0 = interpolant(region.x',region.y'); end