Suggestions for scheme to use to solve PDE numerically

[Chim]

Hello everyone,
I am trying to model the process of laser ablation on a material using MATLAB. The governing equation is of the form:
∂T(x,t)/∂t = ∂/∂x(A*∂T/∂x) + B*exp(-C*t²)*exp(-D*x)
with one Initial condition and two boundary conditions.
Using the built-in 'pdepe' function in Matlab gave inaccurate results, so I have been reading material on solving it numerically. This being a non-linear parabolic pde (correct me if wrong), I am facing difficulty on how to proceed. All the materials that I have referred to so far only address problems where the last term (source term) of the equation is a function of T(x,t). Please guide me.

 

[kai_sikorski]

Is A a function of x or t or T?

If its not a function of any of the above this is just a non-homogenous heat equation and you can solve it analytically.

If it's a function of x or t or both, it's still a linear equation.

If it's a function of T then indeed this is non-linear.

 

[kai_sikorski]

Double Post. Sorry.

 

[bigfooted]

Usually for these problems the Crank-Nicolson method suffices. Note that your source term is quite sensitive to t and x with an exponential decay so a very fine mesh is probably needed when x<<1 and t<<1 to capture the source term accurately.

Also, in numerical methods, it doesn't matter (much) if the source term depends on T or not. In this case, the source term is known on every point of the grid if B,C,D are not functions of T.

You should also be able to get the analytical solution as kai_sikorski mentioned. It is always good to have an analytical solution.

 

[Chim]

Sorry. Forgot to mention about the coefficients earlier.
A is a function of T, whereas B,C and D are positive constants.

Is A a function of x or t or T?

If its not a function of any of the above this is just a non-homogenous heat equation and you can solve it analytically.

If it's a function of x or t or both, it's still a linear equation.

If it's a function of T then indeed this is non-linear.

 

[Chim]

As you have noticed that i will be needing an extremely fine mesh, I was using a time step of 1ns and a 10nm step for distance x. Is it possible that using these increments might affect the stability of the solution?
I would also like to mention that I need to study the behaviour in the first 100ns and the ablation was observed till a few μm.

Usually for these problems the Crank-Nicolson method suffices. Note that your source term is quite sensitive to t and x with an exponential decay so a very fine mesh is probably needed when x<<1 and t<<1 to capture the source term accurately.

Also, in numerical methods, it doesn't matter (much) if the source term depends on T or not. In this case, the source term is known on every point of the grid if B,C,D are not functions of T.

You should also be able to get the analytical solution as kai_sikorski mentioned. It is always good to have an analytical solution.

 

[kai_sikorski]

It might be a good idea to non-dimensionalize time and space. Quantities with units make it hard to tell what's large and what's small. Is 1ns small? Well it depends on how fast the forcing in the problem is varying. If C~(1 ns)^2, or smaller then your time step might be huge in comparison. So I would non-nimensionalize time by (C)^{-1/2}, and space by (D)^{-1}, that way the forcing is changing over order 1 quantities in space and time and the meaning of what a small time step is will be more clear.

Stability is usually related to the time step being too large compared to the space step. I just read the documentation for pdepe though and it seems that it actually chooses an appropriate time step on it's own, the time points you specify are just the time points it will give you a solution at, but they have nothing to do with what time step it actually used. So presumably it chooses something that should be stable.