USER GUIDELINES OF PARABOLIC PDE SOLVER Revision 2

USER GUIDELINES OF: PARABOLIC PDE SOLVER Revision 2. 11 of December 2, 2010 by M. MICCIO

MUC 1. 0 ALLOWS TO SOLVE P. D. E. OF THE KIND : • SECOND ORDER; • PARABOLIC; • LINEAR; • COSTANT COEFFICIENTS.

INITIAL CONDITION : The initial condition can be a general function of x , not only a costant value u 0.

BOUNDARY CONDITIONS : The boundary conditions can be function of time: C(t) e F(t).

WHAT MUC 1. 0 CAN DO: • Resolution of PDE with: explicit Euler method and implicit Crank&Nicholson method. • Check of the Mean Square Error between the solutions.

RESOLUTION OF P. D. E. (1) The parameters to indicate for each specific problem are : Ø Δ; Øk. Ø Lenght of the domain (Lenght); Ø Number of points for discretization (# Points); Ø Time of last solution (Time); Ø Timestep for time discretization (Timestep);

RESOLUTION OF P. D. E. (2) It is necessary to specify initial and boundary conditions in the proper controls. Initial and boundary conditions can be function, respectevely, of x and t.

VALID FUNCTIONS FOR I. C. AND B. C. (1) abs(x) acosh(x) asinh(x) atanh(x) cosh(x) cot(x) Absolute Value Inverse Cosine Inverse Hyperbolic Cosine Inverse Sine Inverse Hyperbolic Sine Inverse Tangent Inverse Hyperbolic Tangent Cosine Hyperbolic Cosine Cotangent

VALID FUNCTIONS FOR I. C. AND B. C. (2) • • csc(x) Cosecant (1/sin(x)) exp(x) Exponential expm 1(x) Exponential (Arg - 1): (e^(x-1)) getexp(x) Mantissa & Exponent (returns the exponent of x) getman(x) Mantissa & Exponent (returns the mantissa of x) int(x) Round To Nearest (rounds its argument to the nearest integer) intrz(x) Round Toward 0 (rounds x to the nearest integer between x and zero) ln(x) Natural Logarithm

VALID FUNCTIONS FOR I. C. AND B. C. (3) • • • lnp 1(x) Natural Logarithm (Arg +1) log(x) Logarithm Base 10 log 2(x) Logarithm Base 2 max(x, y) Maximum min(x, y) Minimum mod(x, y) Quotient & Remainder pow(x, y) x^y rand( ) Random Number (0 - 1) rem(x, y) Remainder sec(x) Secant [computes the secant of x, where x is in radians: (1/cos(x))] • sign(x) Sign (returns 1 if x is greater than 0, returns 0 if x is equal to 0, and returns -1 if x is less than 0)

VALID FUNCTIONS FOR I. C. AND B. C. (4) • sin(x) Sine • sinc(x) Sinc [computes the sine of x divided by x radians: (sin(x)/x)] • sinh(x) Hyperbolic Sin • sqrt(x) Square Root • tan(x) Tangent • tanh(x) Hyperbolic Tangent N. B. 1 The independent variable can be indicated as either x or t N. B. 2 The costant must be indicated as pi(1)

EULER METHOD (1) Timestep and spatial step must be chosen ( ) so that the Stability Parameter be: Otherwise the solution is unstable, and the warning led of the Stability Parameter lights on.

EULER METHOD (2) If the Stability Parameter is greater than 0. 5 the solution is unstable. Pushing the button Euler Method Stability Safety, the s. p. is forced to 0. 499 and the graph does not twist. Anyway the solution plotted is wrong if the Stability Parameter led is red.

CRANK-NICHOLSON METHOD • Implicit Method • Ever stable • Stability Parameter suggested:

Example of students’ test 1) Check the stability of the explicit method for the PDE having = 1 and k = 0 and further subject to a linear initial condition: I. C. : u(x, 0) = 0, 4*x to a Dirichlet condition at the left boundary: B. C. 1: u(0, t) =1 and a a mixed condition at the right boundary: with: Lenght of the domain (Lenght) = 1 Number of points for discretization (# Points) = 40 Time of last solution (Time) = 0. 3 Timestep for time discretization (Timestep) = 0, 00033 2) How much is the spatial step? 3) What is the new value for the time step if the explicit method turns out unstable? 4) Discuss and comment the final diagram 5) Repeat integration with the Crank-Nicholson method and compare the results 15 with the explicit method through the Mean Square Error

STABILITY EXAMPLE 1 – EULER METHOD Parameters: Lenght=1 #Points=40 Time=0, 3 Timestep=0, 00033 =1 k=0 I. C. u(x, o)=0, 4 x A=1 B=0 C=1 D=1 E=-1 F=0 Unstable Solution Initial Condition

STABILITY EXAMPLE 2 – EULER METHOD Parameters: Lenght=1 #Points=40 Time=0, 02 -0, 03 Timestep=0, 00055 =1 k=0 I. C. u(x, 0)=sin(pi(1)*x) A=1 B=0 C=0 D=1 E=0 F=0

EXAMPLE 1 - CRANK & NICHOLSON METHOD Parameters: Lenght=1 #Points=40 Time=1, 7 Timestep=0, 00327024 =1 k=0 I. C. sin(pi(1)*x) A=1 B=0 C=abs(sin(pi(1)*x)) D=1 E=0 F=0

EXAMPLE 2 - CRANK & NICHOLSON METHOD Parameters: L=1 #Points=41 Time=1 Timestep=0, 00082 =1 k=0 I. C. abs(sin(2*pi(1)*x)) A=1 B=0 C=abs(cos(4*pi(1)*x)) D=1 E=0 F=abs(sin(4*pi(1)*x))

COMPARING OF THE SOLUTIONS MUC 1. 0 allows comparison of the solutions calculated with Euler and Crank-Nicholson methods through the function Mean Square Error.

COMPARING OF THE SOLUTIONS (2) The solution vectors – calculated by solving the same problem with each of the two methods – are compared with the MSE appearing in the proper indicator.

MEAN SQUARE ERROR Two vectors having n elements: can be compared calculating the MSE as:

MEAN SQUARE ERROR (2) The MSE found comparing two solutions can not indicate which solution is better: it only calculates a mean error between the solutions. To test the software, some numerical solutions have been compared with analytical solutions or commercial softwares. Here some of the results.

MEAN SQUARE ERROR (3) Parameters: Lenght=1; Time=0, 05; C=1; D=1; E=0; F=1. =1; k=0; I. C. u(x, 0)=0; A=1; B=0; 2 KINDS OF PROOF: • Variable # Points and constant Timestep • Constant # Points and variable Timestep

MUC 1. 0 vs COMMERCIAL SOFTWARE It has been considered the problem of diffusion of a component through a layer of infinite lenght and unitary thickness.

FEM SOLUTION WITH A COMMERCIAL SOFTWARE AUTOMATIC MESH : Number of nodes = 16 Number of elements = 15 Max spatial step = 0, 066667 Time = 0, 1 Timestep = 0, 001

FEM SOLUTION (1)

FEM SOLUTION (2)

MUC 1. 0 SOLUTION Parametri: Lenght=1 # Points=16 Time=0, 1 Timestep=0, 01 =1 k=0 I. C. u(x, 0)=0 A=1 B=0 C=1 D=1 E=0 F=1

MEAN SQUARE ERROR Comparing the Crank & Nicholson Solution Vector with the FEM solution for time 0, 1 : Timestep=0, 01 Timestep=0, 00001 S. P. =1, 1250 S. P. =0, 1125 S. P. =0, 0013 MSE=1, 3 E-3 MSE=3 E-5 MSE= 5, 5 E-6

MUC 1. 0 vs ANALYTICAL SOLUTION (1) Problem of monodimensional diffusion: Analytical solution:

MUC 1. 0 vs ANALYTICAL SOLUTION (2) Lenght=1 # Points=16 Time=0, 1 Timestep=0, 0001 Points X Analytical solution Euler method Crank & Nicholson method 0 0 0 1 0, 667 0, 07753 0, 07777 2 0, 1134 0, 1517 0, 15207 0, 15214 3 0, 2001 0, 2192 0, 21976 0, 21986 4 0, 2668 0, 27708 0, 27784 0, 27798 5 0, 3335 0, 32287 0, 32378 0, 32394 6 0, 4002 0, 35454 0, 35557 0, 35575 7 0, 4669 0, 37069 0, 37182 0, 37200

MUC 1. 0 vs ANALYTICAL SOLUTION (3) MSEA/E = 5, 3 E-7 MSEA/C&N = 7, 3 E-7

Future developments: solving PDE with two and three spacial variables and systems of PDE. PARABOLIC PDE SOLVER by Ugo Avagliano student of Chemical and Food Engineering University of Salerno - Italy