kader <fabdelkader@pi.ac.ae> wrote: > I`m trying to solve the 1D advection-diffusion-reaction equation > dc/dt+u*dc/dx=D*dc2/dx2-kC using Fortan code but I`m still > facing some issues. See: https://class.coursera.org/scientificcomp-003 he does it in Matlab, and it also works in Octave. Actually, that is in 2D, which makes much nicer pictures. You can cheat and go directly to lecture 19, 20, or 21. > first I solved the advection-diffusion equation without > including the source term (reaction) and it works fine. > but when including the source term (decay of substence with > the fisr order decay -kC)I could not get a correct solution. > I have read that I need to use the spliting time mthod or > as called also the operator time method, but I do not know > how to apply it. Are you doing it with Fourier transforms? -- glen

0 |

9/20/2014 3:10:38 AM

On 9/20/2014 5:10 AM, glen herrmannsfeldt wrote: > kader<fabdelkader@pi.ac.ae> wrote: > >> I`m trying to solve the 1D advection-diffusion-reaction equation >> dc/dt+u*dc/dx=D*dc2/dx2-kC using Fortan code but I`m still >> facing some issues. > > See: https://class.coursera.org/scientificcomp-003 Sounds interesting but the link does not seem to work. (The page shows a short *description* of the course but there seems to be no way to actually see it, also not after creating an account.) > he does it in Matlab, and it also works in Octave. > > Actually, that is in 2D, which makes much nicer pictures. > You can cheat and go directly to lecture 19, 20, or 21. No. There is no list of lecture numbers, let alone a way to "goto" any of them. :-) -- Jos

0 |

9/20/2014 8:09:49 AM

Jos Bergervoet <jos.bergervoet@xs4all.nl> wrote: (snip, I wrote) >> See: https://class.coursera.org/scientificcomp-003 > Sounds interesting but the link does not seem to work. > (The page shows a short *description* of the course but > there seems to be no way to actually see it, also not > after creating an account.) I will have to see about that. I can see it because I previously signed up. I am not sure about the schedule. Some will let you sign up anytime, but others won't. -- glen

0 |

9/20/2014 8:32:00 AM

Jos Bergervoet <jos.bergervoet@xs4all.nl> wrote: (snip, I wrote) >> See: https://class.coursera.org/scientificcomp-003 Seems to be: https://class.coursera.org/scientificcomp-006 I didn't realize that the links go to a specific session of the course. The -006 is the current session, that started Aug 26th. Unlike some courses, for this one all homework and quizzes are due on the last day (perfect for procrastinators) which seems to be Oct 31, 11:59PM PDT. (If you set the right time zone, the assignments will do it right. In any case PDT is UTC-0700. You will find that the instructor doesn't follow the discussion forum all that closely, but if you already know about Advection-Diffusion that might not be a problem. You can decide from your background how closely, or not, to follow the lecture, quiz, and homework schedule. There is no cerificate, so it is only for personal gratification and learning. For one of the problems, you do a 2D advection-diffusion equation and watch it move in color. Much of the lecture discussion is on stability of solutions to differential equations. Good stuff to know. -- glen

0 |

9/20/2014 8:48:18 AM

On 9/20/2014 10:48 AM, glen herrmannsfeldt wrote: > Jos Bergervoet<jos.bergervoet@xs4all.nl> wrote: > > (snip, I wrote) > >>> See: https://class.coursera.org/scientificcomp-003 > > Seems to be: https://class.coursera.org/scientificcomp-006 That page contains a text titled "Announcements", which might be a correct title, but there is no actual entry link to the course. From the links in the left-hand pane a tried a few: https://class.coursera.org/scientificcomp-006/wiki/WelcometotheCourse gives again only some description of the course. But then: https://class.coursera.org/scientificcomp-006/lecture finally shows some entry points, only it seems to be unrelated to the numbering you used (in your advice to go to lecture 19, 20, or 21.) Looking, for instance, for the text "Advection-diffusion", I find some sections that might be appropriate, like this one: https://class.coursera.org/scientificcomp-006/lecture/45 I'm not sure whether it is the lecture you actually meant, but it clearly is a video lecture. (I'll leave it to those readers here whit lots of time to further inspect the contents!) .. > Much of the lecture discussion is on stability of solutions > to differential equations. Good stuff to know. It seems to contain a lot about the equations. Does he also treat the solvers? (Stuff like all those different conjugate-gradient implementations, pre-conditioners, multi-grid approach, etc. I don't think I see it listed on the pages I scanned through..) -- Jos

0 |

9/20/2014 9:24:59 AM

Jos Bergervoet <jos.bergervoet@xs4all.nl> wrote: >> Seems to be: https://class.coursera.org/scientificcomp-006 > That page contains a text titled "Announcements", which > might be a correct title, but there is no actual entry > link to the course. From the links in the left-hand > pane a tried a few: > https://class.coursera.org/scientificcomp-006/wiki/WelcometotheCourse > gives again only some description of the course. But then: > https://class.coursera.org/scientificcomp-006/lecture > finally shows some entry points, only it seems to be > unrelated to the numbering you used (in your advice to > go to lecture 19, 20, or 21.) Looking, for instance, for > the text "Advection-diffusion", I find some sections that > might be appropriate, like this one: It seems that the lectures get added week by week. But the lecture notes are available: http://spark-public.s3.amazonaws.com/scientificcomp/Lecture%20Note%20Packet/LectureNotePacketSciComputing.pdf (Maybe only if you are signed up). It says that the lecture notes are for classroom and personal research use only. (snip) > It seems to contain a lot about the equations. Does he > also treat the solvers? (Stuff like all those different > conjugate-gradient implementations, pre-conditioners, > multi-grid approach, etc. I don't think I see it listed > on the pages I scanned through..) It is close to the beginning of solvers, but yes. Well, for Advection-Diffusion you have to solve both the time dependence and space dependence. For time dependence (initial value) he goes through Euler, Runga-Kutta, and Adams, and mostly uses the RK that is part of Matlab and Octave. It is pretty easy, you write the function and pass it to the appropriate routine to call. The spatial solution is done with Fourier transforms. Then to schemes like Leap-Frog (2,2) for combining the time step and spatial solution for a stable result. As always, bigger time step for faster solution, smaller for more accurate solution, but you have to properly consider the time step and spatial mesh size. But much of the discussion is in the direction that the examples and problems are going to use. Oh, the way the assignments work is that you write and run the program, and generate the values of certain time/space points. You then put those in, and, if close enough, you get the points. If you use Octave instead of Matlab, sometimes the nubmers are slightly different, but usually close enough. I would solve it in Octave, then go through and translate to Fortran. The debugging will be a lot faster. -- glen

0 |

9/20/2014 10:12:04 AM

Although we didn't include a source term, which appears to be your primary = difficulty, Chapter 12 of my book has a parallel Burgers equation solver th= at could be a good starting point for you. See the Resources tab at http://= www.cambridge.org/rouson. If you email me, I would be glad to answer quest= ions related to the Burgers equation solver. The Cray, Intel, and GNU Fortran compilers are all capable of compiling it.= (With GNU, you'll need the GCC 5.0, the current development branch. If yo= u want to do parallel runs, you'll also need OpenCoarrays from http://openc= oarrays.org. But you can do serial runs without OpenCoarrays.) Some aspects of a related heat equation solver are described in the modern = Fortran videos at http://jolts.stanford.edu/math.=20 Damian

0 |

9/20/2014 2:43:05 PM

On 9/20/2014 12:12 PM, glen herrmannsfeldt wrote: > Jos Bergervoet<jos.bergervoet@xs4all.nl> wrote: ... >> unrelated to the numbering you used (in your advice to >> go to lecture 19, 20, or 21.) Looking, for instance, for >> the text "Advection-diffusion", I find some sections that >> might be appropriate, like this one: > > It seems that the lectures get added week by week. > > But the lecture notes are available: > > http://spark-public.s3.amazonaws.com/scientificcomp/Lecture%20Note%20Packet/LectureNotePacketSciComputing.pdf Yes, I would say that's the one to start with! > It says that the lecture notes are for classroom and personal research > use only. That's exactly what we need then, it seems. (It saves you hundreds of hours of watching tedious videos! You can clearly see here why the written word is the main foundation of our civilization.) >> It seems to contain a lot about the equations. Does he >> also treat the solvers? (Stuff like all those different >> conjugate-gradient implementations, pre-conditioners, >> multi-grid approach, etc. I don't think I see it listed >> on the pages I scanned through..) > > It is close to the beginning of solvers, but yes. Section 2.3 at least mentions them. But for a real treatment I think "Conjugate Gradient Method Without the Agonizing Pain" still rules supremely: http://www.cs.cmu.edu/~quake-papers/painless-conjugate-gradient.pdf > Well, for Advection-Diffusion you have to solve both the time > dependence and space dependence. Actually Navier-Stokes would be the most interesting example of such a problem (of Millennium-Prize fame! Nathan should give it a dedicated section, I'd say..) > For time dependence (initial value) he goes through Euler, Runga-Kutta, > and Adams, and mostly uses the RK that is part of Matlab and Octave. Which is useful in giving examples of how those Octave routines are called. Saves you hours of reading through tedious help files. > It is pretty easy, you write the function and pass it to the appropriate > routine to call. I'm not sure, however, if we are really coming closer to answering OP's question. (But then again, I forgot what it was..) > ... > I would solve it in Octave, then go through and translate > to Fortran. The debugging will be a lot faster. Well, debugging is only needed if you make errors! So with today's structured programming techniques this hardly seems a factor of importance.. -- Jos

0 |

9/20/2014 4:28:52 PM

On Sunday, September 21, 2014 2:28:57 AM UTC+10, Jos Bergervoet wrote: > On 9/20/2014 12:12 PM, glen herrmannsfeldt wrote: > > > I would solve it in Octave, then go through and translate > > to Fortran. The debugging will be a lot faster. > > Well, debugging is only needed if you make errors! So > with today's structured programming techniques this > hardly seems a factor of importance.. Debugging is always needed, particularly in Fortran, for errors can creep in in many ways. At the very least, runs with test data are needed in order to check that the correct results are obtained. Then there are compiler errors, that may need extensive work to track down and identify. We used to have a customer who, whenever he discovered a programming error in his program, always proclaimed: "Can't understand it! Never made a mistake before in my life!"

0 |

9/21/2014 9:24:48 AM

Op zaterdag 20 september 2014 04:08:34 UTC+2 schreef kader: > Could any one help to provide a simple code for solving this advection-diffusion-reaction equation ? I would appriciate your help. > You have not described the difficulties you experience, but have you thought of reducing the timestep? If the decay coefficient is too large, it will lead to instabilities (i.e. the concentration will grow in size with each time). Alternatively, you can try an implicit method: (C(t+1) - C(t))/Dt + udC(t+1)/dx = Dd2C(t+1)/dx2 - kC(t+1) evaluating all the terms at the new time level. Then the method is stable, although you have to solve a set of linear equations now, but that is straightforward enough. Regards, Arjen

0 |

9/22/2014 9:23:14 AM

Arjen Markus <arjen.markus895@gmail.com> wrote: (snip) > You have not described the difficulties you experience, but > have you thought of reducing the timestep? If the decay > coefficient is too large, it will lead to instabilities > (i.e. the concentration will grow in size with each time). The previously mentioned Coursera course describes the interaction between time step, space step, and instability. > Alternatively, you can try an implicit method: > (C(t+1) - C(t))/Dt + udC(t+1)/dx = Dd2C(t+1)/dx2 - kC(t+1) > evaluating all the terms at the new time level. > Then the method is stable, although you have to solve a set > of linear equations now, but that is straightforward enough. The lectures use a Fourier transform based method for the spatial part. Works well in 2D or even 3D. -- glen

0 |

9/22/2014 9:15:43 PM

Op maandag 22 september 2014 23:15:52 UTC+2 schreef glen herrmannsfeldt: > > The previously mentioned Coursera course describes the > interaction between time step, space step, and instability. > > > > Alternatively, you can try an implicit method:> > > > (C(t+1) - C(t))/Dt + udC(t+1)/dx = Dd2C(t+1)/dx2 - kC(t+1) > > > evaluating all the terms at the new time level. > > Then the method is stable, although you have to solve a set > > of linear equations now, but that is straightforward enough. > > The lectures use a Fourier transform based method for the > spatial part. Works well in 2D or even 3D. > Interesting - I am not familiar with this type of methods. (In most cases the equations I work with have some non-linear reaction term.) I should have a look at this course :). Thanks. Regards, Arjen

0 |

9/23/2014 8:02:14 AM

Arjen Markus <arjen.markus895@gmail.com> wrote: (snip, I wrote) >> The lectures use a Fourier transform based method for the >> spatial part. Works well in 2D or even 3D. > Interesting - I am not familiar with this type of methods. > (In most cases the equations I work with have some > non-linear reaction term.) I should have a look at this > course :). Thanks. http://www.thefouriertransform.com/applications/differentialequations.php (Except much easier to use omega instead of 2*pi*f.) Fourier transforms convert a differential equation into an algebraic equation. Solve that, inverse transform, and you have the solution. The example above, but with omega (w): y''(t)-y(t)=-g(t) Fourier transform both sides, and remember that the Fourier transform of the derivative of a function is -iw times the transform of the function. F{y'}=iw F{y} = iwY So: -w**2*Y(w)-Y(w) = -G(w) solve for Y(w): Y(w)=-G(w)/(1+w**2). Especially convenient for boundary value problems where the function goes to zero at the boundary, or for periodic boundary conditions. -- glen

0 |

9/23/2014 8:34:43 AM

Op dinsdag 23 september 2014 10:34:50 UTC+2 schreef glen herrmannsfeldt: > > Especially convenient for boundary value problems where the function > goes to zero at the boundary, or for periodic boundary conditions. > > -- glen Thanks for this clarification. It has been a while since I last used Fourier transforms or Laplace transforms in this way. Also, I was thinking of Fourier series instead. Regards, Arjen

0 |

9/23/2014 10:38:24 AM

"kader" wrote > Actually I`m using finite difference method to solve the equation and > using the Crank Nicolson scheme for the distretisation. > I`m puting here the methodology of the code but not sure if is it 100 > percent correct specially for the time splitting approach > first is to solve the Advection diffusion equation only and get the system > equation thatcan be solved using TDMA method. > !.....SOLVE EQUATION SYSTEM for the time step DT1 to find the > concentration at time 1/2 > > Ap*C1(i,n+1/2)+AE*C1(i+1,n+1/2)+AW*C1(i-1,n+1/2)=Ap*C1(i,n)+AE*C1(i+1,n)+AW*C1(i-1,n) > CALL TDMA By TDMA do you mean "TriDiagonal Matrix Algorithm"? I recently looked at a different problem in another Fortran related newsgroup that was improved by replacing Gaussian Elimination to solve a linear system by the TDMA method. TDMA takes arguments sub diagonal main diagonal super diagonal right hand side solution vector [size of vectors] diagonals are indexed from row 1, even if it is missing. -- e

0 |

9/26/2014 12:16:20 PM