NR method ODE taylor series R-K method

syms f(x) syms x y(x) f = @(x,y) y^2 - x;

f(x) = input('Enter your function: '); y1 = x - y^2; x = input('Enter the initial value of x: ');
df = diff(f); x0 = 0; y = input('Enter the initial value of y(x): ');
e = input('Enter tolerance: '); y0 = 1; h = input('Enter the step size(h): ');
x0 = input('Enter initial value: '); X = 0.2; X = input('Enter X at which Y is required: ');
n = input('Enter number of iterations: '); h = 0.1; for i = x + h : h : X
if df(x0) ~= 0 y2 = diff(y1); K1 = h * f(x, y);
for i = 1:n y3 = diff(y2); K2 = h * f(x + h/2, y + K1/2);
x1 = x0 - f(x0)/df(x0); y4 = diff(y3); K3 = h * f(x + h/2, y + K2/2);
fprintf('x%d = %.4f\n', i, x1); f1 = diff(y,x,1); K4 = h * f(x + h, y + K3);
if abs(x1 - x0) < e f2 = diff(y,x,2); K = (K1 + 2*K2 + 2*K3 + K4) / 6;
break f3 = diff(y,x,3); x = x + h;
end y10 = subs(y1,{x,y(x)},[x0,y0]); y = y + K;
if df(x1) == 0 y20 = subs(y2,{x,y(x),f1},[x0,y0,y10]); fprintf('Value of y at x = %.1f is %.4f\n', x, y);
disp('N-R method failed'); y30 = subs(y3,{x,y(x),f1,f2},[x0,y0,y10,y20]); end
end y40 = subs(y4,{x,y(x),f1,f2,f3},[x0,y0,y10,y20,y30]); algo: Step 1: Start Step 2: Define function as f(x) Step 3: Define the initial conditions Step 4:
x0 = x1; for i = x0+h:h:X Define the step size Step 5: Compute k1,k2,k3,k4 • k1 is the increment based on the slope at the
end y = y0 + (x - x0)*y10 + ((x - x0)^2/2)*y20 + ((x - x0)^3/6)*y30 + ((x - x0)^4/24)*y40; beginning of the interval, using y • k2 is the increment based on the slope at the midpoint of the
fprintf('\nRoot is: %.4f\n', x1); Y = subs(y,x,i); interval, using y + hk1/2. • k3 is again the increment based on the slope at the midpoint, using y +
else fprintf('Value of y at x = %.1f is %.4f\n', i, Y); hk2/2. • k4 is the increment based on the slope at the end of the interval, using y + hk3 Step 6:
disp('N-R method failed'); end Calculate first degree by Runge-Kutta 4th order Step 7: Increment initial value of x with step size
end algo: Step 1: Start Step 2: Define function as f(x) Step 3. Define the initial conditions Step 4: Step 8: Assign initial value of y as obtained from step 5 Step 9: Display value as required answer
algo: Differentiating each term successively Step 5: Substitution of values in differentiated term Step 10: Stop
Step1: Start Step 2:. Define function as f(x) Step 3: Define first derivative of f(x) as d(f) Step 4: obtained from step 4 Step Step 6: Calculating Taylor series formula by substituting values from
Input initial guess (x0), tolerable error (e) and maximum iteration (N) Step 5: Initialize iteration step 5 Step 7: Display value as required answer Step 8: Stop Milnes predictor
counter i = 1 Step 6: If d(x0) = 0 then print "N R Method Failed" and goto (12) otherwise goto (7)
syms x y
Step 7: Calculate x1 = x0 - f(x0) / g(x0) Step 8: Increment iteration counter i = i + 1 Step 9: If i >= N Modified Eulers f = @(x,y) 1 + x * y^2;
then print "Not Convergent" and goto (12) otherwise goto (10) Step 10: If |f(x1)| > e then set x0 = x1
f = @(x,y) x - y^2; x0 = 0;
and goto (6) otherwise goto (11) Step 11: Print root as x1 Step 12:Stop
x0 = 0; x1 = 0.1;
y0 = 1; x2 = 0.2;
Regula Falsi method xn = 0.4; x3 = 0.3;
syms c f(x) h = 0.1; x4 = 0.4;
f(x) = input('Enter your function: '); x = x0:h:xn; h = 0.1;
a = input('Enter left side of your root guess: '); y(1) = y0; y0 = 1;
b = input('Enter right side of your root guess: '); for i = 1:(length(x)-1) y1 = 1.105;
n = input('Enter the number of iterations you want: '); fprintf('The value is y at (x=%.1f) is = %.4f\n', x(i), y(i)); y2 = 1.223;
e = input('Enter your desired tolerance: '); YE = y(i) + h * f(x(i), y(i)); y3 = 1.355;
if f(a)*f(b) < 0 && a < b for n = 1:3 y11 = f(x1, y1);
for i = 1:n YM = y(i) + (h/2) * (f(x(i), y(i)) + f(x(i+1), YE(n))); y21 = f(x2, y2);
c = (a*f(b) - b*f(a)) / (f(b) - f(a)); YE(n+1) = YM; y31 = f(x3, y3);
fprintf('c%d = %.4f\n', i, c); end y4 = y0 + (4 * h / 3) * (2 * y11 - y21 + 2 * y31);
if f(a)*f(c) < 0 y(i+1) = YM; y41 = f(x4, y4);
b = c; end M = y2 + (h / 3) * (y21 + 4 * y31 + y41);
elseif f(b)*f(c) < 0 fprintf('The value is y at (x=%.1f) is = %.4f\n', xn, YM); fprintf('The value is y at (x=%.1f) is = %.4f\n', x4, M);
a = c;
end algo: Step 1: Start Step2: Define function as f(x) Step 3: Define the initial conditions Step 4: algo: Step 1: Start Step 2: Define function as f(x) Step 3: Define the initial conditions x & Y Step
if abs(f(c)) < e Define step size Step 5: Calculate slope using initial values of x & y) Step 6: Find new y Step 7: 4: Define the step size Step 5: Computing the values of y',y”,y”’ Step 6: Substitute values step 5 in
break Increment x with step size Step 8: Calculate new slope Step 9: Take the mean of two slope Step Milne's predictor formula Step 7: Substitute step 6 value in Milne's corrector formula to get value
end 10: Again find y new , and assign y new with y Step 11; Repeat from step 6 till two consecutive y Step 8: Display value as required answer Step 9: Stop
end are equal Step 12: Increase x and repeat from step 5 till step size Step 13: Display x and
fprintf('\nRoot is: %.4f\n', c); corresponding y Step 14: Stop Gradient
else
syms x y z
disp('No root between given bracket'); Trapezoidal F = x^3 * y^2 * z + y^3 * z^2 * x + z^3 * x^2 * y;
end
syms x disp('Result:')
algo:
f(x) = 1 / (1 + x); gradF = [diff(F,x), diff(F,y), diff(F,z)];
Step 1:Start Step 2: Define function as f(x) Step 3: Read initial left and Right guess (a & b)),
a = input('Enter lower limit a: ');
tolerable error (e) and maximum iteration (N) Step 4: Computer function values f(a) and f(b) Step
b = input('Enter upper limit b: '); algo:step1:start step2:read scalar field step3:calculate 1st order partial derivatives wrt(x,y,z)
5 : Check whether the product of f(a) and f(b) is negative or not. If it is positive take got to step 10.
n = input('Enter the number of sub-intervals n: '); step4:display value of 1st order derivative as gradient step5:substitute the given value of (x,y,z) in
If it is negative then goto step Step 6: for roots from 1 to maximum iteration (N) Determine: c =
h = (b - a) / n; gradient step6:display result step7:stop
[a*f(b) – b*f(a)] / (f(b) – f(a)) Step 7: Display the first root C Step 8: Check whether the value of first
i = 1:1:n-1;
root c is greater than e or not. If yes, goto step 9. If no, goto step 10. (*Here the value e is the desired
degree of accuracy, and hence the stopping criteria.*) Step 9: Check whether the product of f(b)
xi = f(a + (i * h)); Divergence & Curl
I = (h / 2) * (f(a) + f(b) + 2 * sum(xi));
and f(a) is negative or not. If it is negative, then assign a = c; If it is positive, assign b = c; Goto step syms x y z
fprintf('The approximation of above integral is: %f', I);
6 Step 10: Display the root as c. Step 11: Display No root between given bracket Step 12: Stop V = [x^3*y^2*z, y^3*z^2*x, z^3*x^2*y];
a(x,y,z) = V(1);
algo: Step 1: Start Step 2: Read the Function f(x) Step 3: Read the Lower and upper limits (a,b)
Newton forward Step 4: Read the number of sub intervals required Step 5: Calculate width if each interval h = (b –
b(x,y,z) = V(2);
c(x,y,z) = V(3);
x = [100 150 200 250 300 350 400]; a)/n Step 6: Calculate sum of respective terms upto n terms (xi) Step 7: Substitute the obtained
dx = diff(a, x);
y = [10.63 13.03 15.04 16.81 18.42 19.90 21.27]; values in the trapezoidal rule formula (h/2)*(f(a)+f(b)+2*sum(xi))) to find the approximate area of
dy = diff(b, y);
n = length(x); the given curve, Step 8: Display the result Step 9: Stop
dz = diff(c, z);
d = zeros(n-1);
disp('The divergence is:');
h = x(2)-x(1); Simpsons 1/3 div = dx + dy + dz;
xf = 218;
syms x dF = subs(div, [x, y, z], [1, -1, 1]);
p = (xf - x(1)) / h;
f(x) = 1 / (1 + x); colmx = diff(a, z) - diff(c, x);
for k = 2:n
a = input('Enter lower limit a: '); colmy = diff(b, x) - diff(a, y);
d(k-1,1) = y(k) - y(k-1);
b = input('Enter upper limit b: '); colmz = diff(c, y) - diff(b, z);
end
n = input('Enter the number of sub-intervals n: '); disp('The curl is:');
for r = 2:n-1
if rem(n,2) == 1 curl = [colmx; colmy; colmz];
for k = 1:n-r
fprintf('\nEnter valid n!!!'); curlF = subs(curl, [x, y, z], [1, -1, 1]);
d(k,r) = d(k+1,r-1) - d(k,r-1);
n = input('\nEnter n as multiple of 2: ');
end
end algo:step1: start step2:read scalar field step3:calculate vector field with vector components
end
h = (b - a) / n; Fx,Fy,Fz step4:display the vector field step5:calculate 1st order partial derivatives wrt (x,y,z)
disp('The difference table is:')
i = 1:1:n-1; step6:calculate divergence as sum of values obtained from step 5 step7:substitute the given
d
X = f(a + i * h); value of (x,y,z) in divergence step8:calculate components of curl vector values obtained from step
s = y(1);
even = sum(X(2:2:n-1)); 5 step9:substitute the given values of(x,y,z) in curl step10:display result step11:stop
t = p;
odd = sum(X(1:2:n-1));
for r = 1:n-1
s = s + t * d(1,r);
I = (h / 3) * (f(a) + f(b) + 4 * odd + 2 * even); Inverse laplace
fprintf('The approximation of above integral is: %f ', I);
t = (p - r)/(r + 1) * t; syms s t
end F1 = 1/s;
algo: Step 1: Start Step 2: Read the Function f(x) Step 3: Read the Lower and upper limits (a,b)
fprintf('The required value is f(%.1f) = %.4f\n', xf, s); F2 = 1/s^2;
Step 4: Read the number of sub intervals required Step 5: Calculate width if each interval h = (b –
F3 = 1/(s^2 + 1);
a)/n Step 6 : Calculate sum of odd and even terms separately from the sub- intervals upto n terms
algo: Step 1: Start Step 2: Read the n equidistant x values and corresponding y values into arrays F1 = ilaplace(F1, s, t);
Step 7: Substitute all these values in Simpson's 1/3 rule formula(b∫ₐ f(x) dx =(h/3)
x and y. Step 3: Calculate the length of array x Step 4: Initialize a matrix “D” to store divided F2 = ilaplace(F2, s, t);
[(f(a)+f(b))+4*Odd)+2*even) and simplify to find the approximate area of the given curve, Step 8:
differences. Step 5: Calculate the common difference, h, between consecutive x values: h = x[1] F3 = ilaplace(F3, s, t);
Display the result Step 9: Stop
- x[0] Step 6: Read the value of x_interp, the point where interpolation is desired. Step 7: Initialize disp(['Inverse Laplace transform of 1/s: ', char(F1)]);
the result to (xf-x(1)/h (the first divided difference). After that, initialize term to 1.0. Step 8: disp(['Inverse Laplace transform of 1/s^2: ', char(F2)]);
Calculation of first forward differences using Newton’s forward difference formula Step 9: Simpsons 3/8 disp(['Inverse Laplace transform of 1/(s^2 + 1): ', char(F3)]);
Calculation of second and rest forward differences Step 12: Display the difference table “D” Step syms x
13:Calculate value of x_interp, the point where interpolation is Step 14: Display X-Interp Value f(x) = 1 / (1 + x.^2); algo: Step 1: Define symbolic variables Step 2: Read Standard Laplace transforms Step 3:
a = input('Enter lower limit a: '); Compute inverse Laplace transforms Step 4: Display the results Step 5: Stop
Newton backward b = input('Enter upper limit b: ');

x = [80 85 90 95 100];
n = input('Enter the number of sub-intervals n: '); Laplace Transform of Convolution
h = (b - a) / n;
y = [5026 5674 6362 7088 7854]; syms s t
xi = a:h:b;
xf = input('Enter the value of x where we want to find the value of f(x): '); f = exp(-t);
yi = f(xi);
n = length(x); g = t^2;
i = 1:1:n-1;
d = zeros(n-1); F = laplace(f, s);
S = f(a + i * h);
h = x(2) - x(1); G = laplace(g, s);
I = 3:3:n-1;
p = (xf - x(n)) / h; H = F * G;
S3 = sum(S(I));
for k = 2:n Convolution = ilaplace(H, s, t);
S(I) = [];
d(k-1,1) = y(k) - y(k-1); disp(['Laplace transform of convolution: ', char(H)]);
SO = sum(S);
end disp(['Inverse Laplace transform of convolution:', char(Convolution)]);
I = (3 * h / 8) * (f(a) + f(b) + 3 * SO + 2 * S3);
for r = 2:n-1
fprintf('The approximation of above integral is %f\n', I);
for k = 1:n-r algo: Step 1: start Step 2: Define symbolic variables Step 3: Read the first and second functions
d(k,r) = d(k+1,r-1) - d(k,r-1); Step 4: Compute inverse Laplace transforms Step 5: Compute the Laplace transform of the
algo: Step 1: Start Step 2: Read the Function f(x) Step 3: Read the Lower and upper limits (a,b)
end convolution Step 6: Compute the inverse Laplace transform of the convolution Step 7: : Display
Step 4: Read the number of sub intervals required Step 5: Calculate width if each interval h = (b –
end the results Step 8: stop
a)/n Step 6: Calculate sum all polynomial with degree 3 terms (S3) and delete the number from
disp('The backward difference table is:')
the array Step 7: Calculate sum of remaining numbers (S0) Step 8: Substitute all these values in
d
Simpson's 3/8 rule formula(b∫ₐ f(x) dx = (3*h/8)*((f(a)+f(b)+3*SO + 2*S3))and simplify to find the
s = y(n);
approximate area of the given curv Step 9: Display the result Step 10: Stop
t = p;
m = n - 1;
for r = 1:n-1
s = s + t * d(m,r);
t = (p + r)/(r + 1) * t;
m = m - 1;
end
fprintf('The required value is f(%.0f) = %.4f\n', xf, s);

algo:
Step 1: Start Step 2: Read the n equidistant x values and corresponding y values into arrays x and
y. Step 3: Calculate the length of array x Step 4: Initialize a matrix “D” to store divided differences.
Step 5: Calculate the common difference, h, between consecutive x values: h = x[1] - x[0] Step 6:
Read the value of x_interp, the point where interpolation is desired. Step 7: Initialize the result to
(xf-x(n)/h (the first divided difference). After that, initialize term to 1.0 Step 8: Calculation of first
backword differences using Newton’s forward difference formula Step 9: Calculation of second
and rest backword differences Step 10: Display the difference table “D” Step 11: Calculate value
of x_interp, the point where interpolation is Step 12: Display X-Interp Value