How to handle interface boundary condition in convection-diffusion reaction equation for pressure, velocity and concentration equation?

Complete mathematical model equations and boundary condition is given in the attached pdf. the domain for the problem taken as 2D disk with three region separated by different radius (given the attached pdf). I am trying to solve the problem using finite element method using MATLAB Partial Differential Equation Toolbox ( MATLAB pdetoolbox). There in the boundary condition section it is given that for the inside regions it follows interface boundary condition where it gives continuity for the variables (see the attached pdf boundary condition). but in the pdetoolbox manual it is written, we can not implement boundary condition for subdomains boundaries. but there must be some way to handle these interface boundary condition, may be by making these equation as coupled equation with the original equation, or some otherway in the specify boundary condition. I want to solve the problem using finite element method, and very new to pdetoolbox. i am seeking help from the community, please help me out.

although i tried to solve the problem using dirichlet and neumann boundary condition but answer and graph not coming as expected. I am providing some portion of my code, please help me out.

G = [0;0;0];

Q=[0;0;0;0;0;0;0;0;0];

u0=[0;0;10];

applyBoundaryCondition(model,"dirichlet","Edge",9:12,"u",[0;0;0]);

applyBoundaryCondition(model,"neumann","Edge",1:4,"g",G,"q",Q);

applyBoundaryCondition(model,"neumann","Edge",5:8,"g",G,"q",Q);

setInitialConditions(model,u0);

specifyCoefficients(model,"m",0,"d",dn,"c",cn,"a",an,"f",fn,"Face",3);

specifyCoefficients(model,"m",0,"d",dt,"c",cv,"a",av,"f",fv,"Face",2);

specifyCoefficients(model,"m",0,"d",dh,"c",ch,"a",ah,"f",fh,"Face",1);

Answer 1

[Torsten]

If the interface conditions are continuity of the variable and the flux, the finite element method will automatically respect them - you don't need to set any conditions at the internal boundaries. If you want other interface conditions (that might perhaps cause discontinuous solutions or something exotic), you cannot use the PDE toolbox.

[Gobinda Debnath]

thank you so much for your answer, and what about the boundary condition that says

how can we implement this?

[Torsten]

What problem do you encounter when you try to specify the above boundary condition ?

Is your geometry cylindrical or spherical ? Did you pay attention to this in the setup of your problem ?

[Gobinda Debnath]

Dear sir, I was giving proper attention to the problem i have read and coded the entire problem and understand it, and i understand the interface boundary condition as it maintained autamatically by pde toolbox, some queries and doubts are as follows

1. as i have asked in last comments, in pde toolbox we can give boundary condition at outer edge and if there is hole then inside boundary we can give boundary condition, but at radius r=o that is at the centre of sphere how can give boundary condition?
2. As in the attched PDF file i have given the problem there you may have seen in equation (2), then how to handle this term in pde toolbox? I have considered it for 2D case as , then used function handle to code this, am I correct?
3. in the attched PDF file i have given the problem there you may also have seen in equation (3), then how to handle this term?, since here divergence operator and two variable, therefore, I have considered it for 2D case as

then I used function handle to tackle this in PDE toolbox, Am I correct on this, Please help me on this 3 doubts thankyou very much.

[Torsten]

Since you are talking about r = 0 several times now, it seems you are working in cylindrical or spherical coordiantes.

In this coordinate system, r = 0 is an outer boundary and you can specify your boundary condition as du/dr = 0, I guess. But notice that the divergence operator in cylindrical or spherical coodinates is different from the cartesian. Thus div(nu*C) as written above looks incorrect.

[Gobinda Debnath]

hello sir, @Torsten please attached pdf see the equation 2, there is v=-K*delP, that is gradient operator now how can we write f coefficient for this variable v.

in the above pdf three variable are there v,P,C, so for first variable v, I have given f coefficient in this way but it shows error,

i just want to know how to handle gradient operation in math equation in matlab pde toolbox, please help

fn=@(location,state)[(-knecrotic*gradient(u(2,:))).*(location.x+location.y);0+0*(location.x+location.y);-(state.u(1,:).*state.ux(3,:)+state.u(1,:).*state.uy(3,:)+state.u(3,:).*state.ux(1,:)+state.u(3,:).*state.uy(1,:))];

fv=@(location,state)[(-ktumor*gradient(u(2,:))).*(location.x+location.y);pss1+0*(location.x+location.y);-(state.u(1,:).*state.ux(3,:)+state.u(1,:).*state.uy(3,:)+state.u(3,:).*state.ux(1,:)+state.u(3,:).*state.uy(1,:))];

fh=@(location,state)[(-khealthy*gradient(u(2,:))).*(location.x+location.y);pss2+0*(location.x+location.y);-(state.u(1,:).*state.ux(3,:)+state.u(1,:).*state.uy(3,:)+state.u(3,:).*state.ux(1,:)+state.u(3,:).*state.uy(1,:))];

[Torsten]

I have not enough experience with the PDE toolbox to answer your question.

But why don't you use "pdepe" for your problem ? It seems to be one-dimensional with only r as independent spatial variable.

[Hau Trien]

Hi Torsten,

In case I use different type of elements on different zones, will continuity interface conditions apply automatically? For example, I use P2-P1 on zone 1 and RT0-P0 on zone 2 for velocity and pressure, respectively.

[Torsten]

The forum has little experience with the application of the PDE Toolbox to Fluid Flow Problems.

Better ask Technical Support for such specific questions:

https://uk.mathworks.com/support/contact_us.html