How can the pdepe solver solution be used as a new initial condition for the next iteration? i.e using pdepe within a loop

조회 수: 13 (최근 30일)
I solved the convection diffusion equation using pdepe solver.
I would like to solve this equation within a loop. i.e
the output solution at the tout is a vector, which should be used as new initial condition.
It would be appreciated to help in this issue.
  댓글 수: 1
maria atiq
maria atiq 2024년 1월 30일
Hey, were you able to solve this issue? I want to use my pdepe solver within a loop where i have designed a PID to control the output of my pde. Kindly share if you were able to run the pdepe as a loop.

댓글을 달려면 로그인하십시오.

답변 (1개)

Venkatachala Sarma
Venkatachala Sarma 2016년 1월 12일
편집: Venkatachala Sarma 2016년 1월 12일
Hi Younus,
This could be done by using a variable shared between the function calling 'pdepe' and the function pdex1ic(x) which gets called to manage the initial conditions. Let us take a global variable as an example.
Consider the first example in the following documentation page.
In this example, you could declare a global variable which will be updated with the solution of 'pdepe' every time it runs. Just equate this global variable to the output of the function pdex1ic which manages the initial condition.
I am pasting the modified version of pdex1 and pdex1ic functions alone. Remaining of the example code is retained.
if true
function pdex1
global takeOut;
m = 0;
x = linspace(0,1,20);
t = linspace(0,2,5);
takeOut = 0;
for i = 1:5
sol = pdepe(m,@pdex1pde,@pdex1ic,@pdex1bc,x,t);
u = sol(:,:,1);
takeOut = u(end,end);
title('Numerical solution computed with 20 mesh points.')
xlabel('Distance x')
ylabel('Time t')
title('Solution at t = 2')
xlabel('Distance x')
The above function puts the pdepe function in a loop and uses the global variable named takeOut to carry the current solution. This is just to show an example and you can explore the other possibilities. The pdex1ic function shown below just copies the value in this variable.
if true
function u0 = pdex1ic(x)
global takeOut;
u0 = takeOut;
Update these functions in the example and run it. You could see that the outputs vary depending upon the initial conditions.
Hope this helps.
Regards Venkatachala Sarma
  댓글 수: 4
Antonio López
Antonio López 2018년 8월 23일
Dear Andreas,
I've found a solution.
Let me denote x as the spatial coordinate. For example, imagine that x is disretised between 0 and 0.05 with intervals of 0.01, so x=0:0.01:0.05. The lenght of this vector is 6.
From the solution of our first PDE, we obtain u1. For example, it could be that u1(end,:)=[1,2,3,4,3,1] (I'm inventing the numbers). Of course, the length of u1(end,:) is 6.
First of all, you have to copy your spatial coordinate discretization to another variable, for example, xcopy. This variable must be declared as global so that it could be accesible from the initial condition function.
In the initial condition function of the next PDE, you can do this:
function u0 = pdexAAAAA_ic(x)
global xcopy;
u0 = u1(end, find(xcopy==x,1,'first'));
The goal of this initial condition function is just to fill the first row of the next PDE solution, so internally it's like a for loop in which x is the index of the loop. Every time u0 is called, you are finding the position of the current value of x in the vector xcopy. Doing this, you can access the previous solution u1 in the time end and in the correct position.
I hope it would be helpful. Inserting a breakpoint in the initial condition function would help to understand the process during the execution.
Haonian Shu
Haonian Shu 2020년 4월 20일
Dear Antonio,
I am have the same issue and have been trying to follow your solution. Yet I cannot figure out where should I declare the solution u1, since u1 is not global, how can it be extracted by the ic function. also, I am struggling to understand the goal of this initail condition function, the goal of xcopy. actually, all i want to find is to give pdepe a custumized, discrete (arbitrary) initial condition rather than a number or a simple analytical function.

댓글을 달려면 로그인하십시오.


Help CenterFile Exchange에서 Boundary Conditions에 대해 자세히 알아보기

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by