Subodh Acharya
2010-Apr-26 01:09 UTC
[R] finite difference scheme for 2D differential equations
Hello everyone, I am trying to solve 2D differential equations using finite difference scheme in R. I have been able to work with the equations with only one spatial dimensions but I want to extend it to the two dimensional problem. For example i can simulate one dimensional diffusion using a code like the following. But I want to write a similar code for,say, a two dimensional diffusion equation. Any kind of help/advice is highly appreciated. Here is how I did for the 1D equation (in this case for 1D diffusion) Equation: dc/dt = D*d^2c/dx^2 dx = 10 dt = 0.1 D = 15 n = 20 tstep = 200 C = 50 dt = 0.1 dx = 10 conc<-matrix(0,tstep,n) conc[1,1:10] = C for (i in 2:tstep) { for (j in 1:n){ if (j==1){ conc[i,j] = conc[i-1,j] - dt*(D/dx^2)*(conc[i-1,j] - conc[i-1,j+1]) } conc[i,j] = conc[i,j] if(j >1 && j < n){ conc[i,j] = conc[i-1,j] + dt*(D/dx^2)* (conc[i-1,(j-1)] - 2*conc[i-1,j] + conc[i-1,(j+1)]) } if (j==n){ conc[i,n] = conc[i-1,n] + dt*(D/dx^2)* ( conc[i-1,n-1] - conc[i-1,n]) } conc[i,j] = conc[i,j] } } Now in 2D the equation will be like this dc/dt = Dx*d^2c/dx^2 + Dy*d^C/dy^2 So that when I solve it I will get C(x, y, t) at each node of the grid. Thanks for the help in advance. -- Acharya, Subodh [[alternative HTML version deleted]]
Thomas Petzoldt
2010-May-12 06:25 UTC
[R] finite difference scheme for 2D differential equations
Hello Subodh Acharya, I've been away for a field trip for two weeks, and I guess you have already found a solution for your problem. If your question is still open, you may consider to look into package deSolve that now provides functions (ode1D, ode2D and ode3D) for solving 1-3 dimensional equations using the method of lines approach. The package comes with many examples and extensive documentation and there is also a JSS publication: http://www.jstatsoft.org/v33/i09 If you need more help, please write to the dedicated mailing list for dynamic models: https://stat.ethz.ch/mailman/listinfo/r-sig-dynamic-models Hope it helps Thomas Petzoldt PS: You may also have a look into package ReacTran. At 26.04.2010 03:09 Subodh Acharya wrote:> Hello everyone, > I am trying to solve 2D differential equations using finite difference > scheme in R. I have been able to work with the equations with only one > spatial dimensions but I want to extend it to the two dimensional problem. > For example i can simulate one dimensional diffusion using a code like the > following. But I want to write a similar code for,say, a two dimensional > diffusion equation. Any kind of help/advice is highly appreciated. > > Here is how I did for the 1D equation (in this case for 1D diffusion) > Equation: dc/dt = D*d^2c/dx^2 > > > dx = 10 > dt = 0.1 > D = 15 > n = 20 > tstep = 200 > C = 50 > dt = 0.1 > dx = 10 > conc<-matrix(0,tstep,n) > conc[1,1:10] = C > for (i in 2:tstep) { > for (j in 1:n){ > if (j==1){ > conc[i,j] = conc[i-1,j] - dt*(D/dx^2)*(conc[i-1,j] - conc[i-1,j+1]) > } > conc[i,j] = conc[i,j] > if(j>1&& j< n){ > conc[i,j] = conc[i-1,j] + dt*(D/dx^2)* (conc[i-1,(j-1)] - 2*conc[i-1,j] + > conc[i-1,(j+1)]) > } > if (j==n){ > conc[i,n] = conc[i-1,n] + dt*(D/dx^2)* ( conc[i-1,n-1] - conc[i-1,n]) > } > conc[i,j] = conc[i,j] > } > } > > Now in 2D the equation will be like this > > dc/dt = Dx*d^2c/dx^2 + Dy*d^C/dy^2 > > So that when I solve it I will get C(x, y, t) at each node of the grid. > > Thanks for the help in advance.