LUtri modul za rjesavanje trodijagonalnog sustava linearnih jednadzbi

LUtri[n_, a_, b_, c_, y_] := 
      Module[{d, e, f, z, x}, 
 (* 
          Treba unijeti red matrice n,  vektore glavne a,  
          gornje b i donje c dijagonale,  
          te vektor slobodnih koeficijenata y *) 
 
        d =   e =  f = z =  x = 
                  Table[0, {i, 
                      n}]; 
                                            (* 
          Izracunavanje dekompozicije *) 
e[[1]] = 
          a[[1]]; f[[1]] = b[[1]]; 
        d[[1]] = c[[1]]/e[[1]];  
        Do[e[[i]] = a[[i]] - d[[i - 1]] f[[i - 1]]; 
          f[[i]] = b[[i]];  
          d[[i]] = c[[i]]/e[[i]], {i, 2, 
            n - 1}]; 
e[[n]] = 
          a[[n]] - 
            d[[n - 1]] f[[n - 
                    1]]; 
                                            (* 
          Rjesavanje sustava LUx = y *) 
z[[1]] = 
          y[[1]]; 
Do[
          z[[i]] = y[[i]] - d[[i - 1]] z[[i - 1]], {i, 2, 
            n}]; 
x[[n]] = 
          z[[n]]/e[[n]]; 
 
        Do[x[[n - i]] = (z[[n - i]] - 
                  f[[n - i]] x[[n - i + 1]])/e[[n - i]], {i, 
            n - 1}]; 
x // N]

Primjer 1:  y''=-y;     y[0]=3,    y[Pi/2]=7

Clear[x, y]
DSolve[{y''[x] == -y[x], y[0] == 3, 
        y[Pi/2] == 7}, y[x], x]

yy[t_] := 3  Cos[t] + 7  Sin[t];
A = 0;  B = Pi/2;  alpha = 3;  beta = 7;
u[x_] := 0;  v[x_] := -1;  w[x_] := 0;
n = 
      10; h = (B - A)/(n + 
              1) ;

Primjer 2: y''=2Exp[x]-y;      y[0]=2,    y[1]=E+Cos[1],    sol[t]=3Cos[t]+7Sin[t]

Clear[x, y]
DSolve[{y''[x] == 2  Exp[x] - y[x], y[0] == 2, 
        y[1] == E + Cos[1]}, y[x], x]

yy[x_] := 
      Cos[x] + (E^x) Cos[x]^2 - E Cos[1] Cot[
            1] Sin[x] + E Csc[1] Sin[
            x] - E Sin[1] Sin[
            x] + (E^x) Sin[x]^2
A = 0;  B = 1;  alpha = 2;  
    beta = E + Cos[1];
u[x_] := 2  Exp[x];  v[x_] := -1;  
    w[x_] := 0;
n = 
      10; h = (B - A)/(n + 
              1) ;

Metoda konacnih diferencija

x = Table[A + i h, {i, 0, n + 1}];
a = 
        Table[2 + h^2  v[x[[i]]], {i, n}];
b = 
        Table[-1 +  .5  h w[x[[i]]], {i, 
            n}];
c = 
        Table[-1 -  .5  h w[x[[i + 1]]], {i, 
            n}];
d = Table[-h^2 u[x[[i]]], {i, n}];  
    d[[1]] = d[[1]] + (1 +  .5  h w[x[[1]]]) alpha;  
    d[[n]] = d[[n]] - b[[n]] beta;
y = LUtri[n, a, b, c, d];
slf = 
        Plot[yy[t], {t, A, B}, PlotRange -> {0, 8}, 
          DisplayFunction -> Identity];
rub = 
        ListPlot[{{A, alpha}, {B, beta}}, 
          PlotStyle -> {Hue[1], AbsolutePointSize[5]}, 
          DisplayFunction -> Identity];
sol = 
        ListPlot[Table[{x[[i + 1]], y[[i]]}, {i, n}], 
          PlotStyle -> {Hue[ .7], AbsolutePointSize[4]}, 
          DisplayFunction -> Identity];
Show[slf, rub, sol, 
        DisplayFunction -> $DisplayFunction];

Primjena gotovih Mathematica modula

Clear[x, y]
DSolve[{y''[x] == -y[x], y[0] == 3, 
        y[Pi/2] == 7}, y[x], x]
Plot[Evaluate[y[x] /. %], {x, 0, Pi/2}, 
        PlotRange -> {0, 8}];