Solve Bratu problem using Python

784 Views Asked by At

I am going crazy trying to solve the Bratu problem using Python:

$$y''(x)+ e^{y(x)} = 0, \quad \lambda = 1, \quad x \in(0,1),$$

$$y(0) = y(1) = 0$$

I have to solve this using the tridiagonal matrices. None of the solutions that I have written work so far. Any help in this direction??

1

There are 1 best solutions below

1
On

This seems to get to (one of them) the solution:

clear all

% Initialize the solution

yl      = 0;
yr      = 0;
lambda  = 1;

% Set up the parameters:

N       = 200;
x0      = 0;
xf      = 1;
dx      = (xf-x0)/N;
x       = x0:dx:xf;

y       = zeros(N+1,1);
y(1)    = yl;
y(end)  = yr;
dY      = zeros(N+1,1);

tol     = 1e-12;
kmax    = 500;
err     = 1;
k       = 0;

C       = zeros(N-2,1);

% Begin NR iteration:

while(err>tol && k<kmax)

    d   = -(2+dx^2*lambda*exp(y(2:N)));
    infd= ones(length(d)-1,1);
    supd= ones(length(d)-1,1);

    A   = tridiagonal(d,infd,supd);
    for i = 2:N
       C(i-1) = -( y(i+1)-2*y(i)+y(i-1) + dx^2*lambda*exp(y(i))); 
    end

    dY(2:N)  = A\C;
    y        = dY + y;    

    k       = k+1;
    err     = norm(dY,2)/norm(y,2);   

    fprintf('Error between iterations: %f \n',err);    

    plot(x,y,x,dY,'-.');    
    xlabel('x');
    ylabel('y');
    drawnow;

end

where the function 'tridiagonal' is given by:

function A = tridiagonal(d,infd,supd)
A = diag(d) + diag(infd,-1) + diag(supd,1);
end

This is a MATLAB code using a usual Newton-Raphson method to approach the solution. If you are interested only in coding, I assume you are familiar enough with this and can easily translate my code to Python (which it's almost straightforward).

If you are also interested in the mathematical steps behind this, please let me know.

Hope this helps.

Cheers!