solving a second order non linear differential equation using RK 4TH order method

2023 4 3
2 个回答
更新时间：2025 2 5
4 次查看（30 天）

0 个投票

Differential equation : h d^2h/dx^2 + (dh/dx)^2 - dh/dx * tan(ax) + c - h * sec^2(ax) * a = 0
Boundary conditions: h(x=0)=h0 and h(x=L)=h0
Dependent variable: h
Independent variable: x
constants: a,c,L,h0
Method to be used : RK 4th order
please help me
let y = h and z = dh/dx=dy/dx,dz/dx = a * sec^2(ax) + (1/h) * (z * tan(ax) - z^2 - c) confused how to give boundary conditions

6 个评论
显示 4更早的评论 隐藏 4更早的评论

Jack
Jack 2023-4-3
To use the fourth-order Runge-Kutta method to solve this differential equation with the given boundary conditions, you can use the following steps:
  1. Define the dependent variable y as h and the independent variable x as x.
  2. Define the initial conditions y(0) = h0 and z(0) = 0.
  3. Define the constants a, c, and L.
  4. Define the step size dx as a small positive value.
  5. Set up a loop to iterate over x values from 0 to L, with a step size of dx.
  6. Inside the loop, use the current values of y and z to calculate the derivatives dy/dx and dz/dx using the given differential equation.
  7. Use the fourth-order Runge-Kutta method to update the values of y and z at the next x value.
  8. Check if the current x value is equal to either 0 or L. If so, update the value of y at that boundary condition.
  9. Repeat the loop until the final x value is reached.
Here is the code that implements the above steps:
function [x, y] = solve_differential_equation(a, c, L, h0)
% Solve the differential equation h''(x) + (h'(x))^2 - h'(x) * tan(a*x) + c - h(x) * sec^2(a*x) * a = 0
% with boundary conditions h(0) = h0 and h(L) = h0, using the fourth-order Runge-Kutta method.
% Define initial conditions
y0 = h0; % h(0) = h0
z0 = 0; % h'(0) = 0
% Define step size
dx = 0.01;
% Set up loop
x = 0:dx:L;
n = length(x);
y = zeros(n, 1);
z = zeros(n, 1);
y(1) = y0;
z(1) = z0;
for i = 1:n-1
% Calculate derivatives
dydx = z(i);
dzdx = a * sec(a*x(i))^2 + (1/y(i)) * (z(i) * tan(a*x(i)) - z(i)^2 - c);
% Update y and z using fourth-order Runge-Kutta method
k1 = dx * dydx;
l1 = dx * dzdx;
k2 = dx * (dydx + 0.5*l1);
l2 = dx * (dzdx + 0.5*k1);
k3 = dx * (dydx + 0.5*l2);
l3 = dx * (dzdx + 0.5*k2);
k4 = dx * (dydx + l3);
l4 = dx * (dzdx + k3);
y(i+1) = y(i) + (1/6)*(k1 + 2*k2 + 2*k3 + k4);
z(i+1) = z(i) + (1/6)*(l1 + 2*l2 + 2*l3 + l4);
% Check boundary conditions
if x(i+1) == 0
y(i+1) = h0;
elseif x(i+1) == L
y(i+1) = h0;
end
end
% Plot results
plot(x, y)
xlabel('x')
ylabel('h(x)')
title('Solution to differential equation')
end
To use this function, call it with the values of a, c,
Torsten
Torsten 2023-4-3
RK4 is a solution method usually used for initial value problems. Yours is a boundary value problem. How do you want to use RK4 ? Together with single or multiple shooting ?
Torsten
Torsten 2023-4-6
编辑：Torsten 2023-4-6
And why don't you use MATLAB's bvp4c instead of progamming a complicated coupling of a root finder and an integrator on your own ? Masochism ?
Can you please help me with bvp4c or bvp5c i am kind of stuck Thank you

请先登录，再进行评论。

请先登录，再回答此问题。

回答（2 个）

Jack
Jack 2023-4-3

0 个投票

Here's an example code in MATLAB for solving the given differential equation using the RK4 method. Note that you need to define the constants and initial conditions before running the code.
% Define constants and initial conditions
a = 1;
c = 1;
L = 1;
h0 = 1;
N = 1000; % Number of grid points
x = linspace(0, L, N)';
dx = x(2) - x(1);
h = h0*ones(N, 1); % Initial guess for h
dhdx = zeros(N, 1); % Initial guess for dh/dx
% Define the function f(x, y) = [dy/dx, d^2y/dx^2]
f = @(x, y) [y(2); ...
(-y(2)^2 + y(2)*tan(a*x) - c + h0*sec(a*x)^2*a)/h0];
% Solve the differential equation using the RK4 method
for n = 1:10000
k1 = dx*f(x, [h, dhdx]);
k2 = dx*f(x + dx/2, [h + k1(1:N)/2, dhdx + k1(N+1:end)/2]);
k3 = dx*f(x + dx/2, [h + k2(1:N)/2, dhdx + k2(N+1:end)/2]);
k4 = dx*f(x + dx, [h + k3(1:N), dhdx + k3(N+1:end)]);
h = h + (k1(1:N) + 2*k2(1:N) + 2*k3(1:N) + k4(1:N))/6;
dhdx = dhdx + (k1(N+1:end) + 2*k2(N+1:end) + 2*k3(N+1:end) + k4(N+1:end))/6;
end
% Plot the solution
plot(x, h);
xlabel('x');
ylabel('h');
title('Solution of the differential equation');

类别

在 帮助中心 和 File Exchange 中查找有关 Numerical Integration and Differential Equations 的更多信息

标签

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by