Newton's Method on MATLAB for Stationary Solutions for the Non-linear Klein-Gordon Equation
8 views (last 30 days)
Show older comments
Athanasios Paraskevopoulos
on 17 Mar 2024
By interpreting the equation in this way, we can relate the dynamics described by the discrete Klein - Gordon equation to the behavior of DNA molecules within a biological system .
I am trying to represnt stationary points of the discrete Klein - Gordon equation. And the result are as follows.
Mathematica are different from MATLAB Any suggestion?
UPDATED CODE
% Parameters
numBases = 100; % Number of spatial points
omegaD = 0.2; % Common parameter for the equation
% Preallocate the array for the function handles
equations = cell(numBases, 1);
% Initial guess for the solution
initialGuess = 0.01 * ones(numBases, 1);
% Parameter sets for kappa and beta
paramSets = [0.1, 0.05; 0.5, 0.05; 0.1, 0.2];
% Prepare figure for subplot
figure;
set(gcf, 'Position', [100, 100, 1200, 400]); % Set figure size
% Newton-Raphson method parameters
maxIterations = 1000;
tolerance = 1e-10;
for i = 1:size(paramSets, 1)
kappa = paramSets(i, 1);
beta = paramSets(i, 2);
% Define the equations using a function
for n = 2:numBases-1
equations{n} = @(x) -kappa * (x(n+1) - 2 * x(n) + x(n-1)) - omegaD^2 * (x(n) - beta * x(n)^3);
end
% Boundary conditions with specified fixed values
someFixedValue1 = 10; % Replace with actual value if needed
someFixedValue2 = 10; % Replace with actual value if needed
equations{1} = @(x) x(1) - someFixedValue1;
equations{numBases} = @(x) x(numBases) - someFixedValue2;
% Combine all equations into a single function
F = @(x) cell2mat(cellfun(@(f) f(x), equations, 'UniformOutput', false));
% Newton-Raphson iteration
x_solution = initialGuess;
for iter = 1:maxIterations
% Evaluate function
Fx = F(x_solution);
% Evaluate Jacobian
J = zeros(numBases, numBases);
epsilon = 1e-6; % For numerical differentiation
for j = 1:numBases
x_temp = x_solution;
x_temp(j) = x_temp(j) + epsilon;
Fx_temp = F(x_temp);
J(:, j) = (Fx_temp - Fx) / epsilon;
end
% Update solution
delta = -J\Fx;
x_solution = x_solution + delta;
% Check convergence
if norm(delta, inf) < tolerance
fprintf('Convergence achieved after %d iterations.\n', iter);
break;
end
end
if iter == maxIterations
error('Maximum iterations reached without convergence.');
end
% Plot the solution in a subplot
subplot(1, 3, i);
plot(x_solution, 'o-', 'LineWidth', 2);
grid on;
xlabel('n', 'FontSize', 12);
ylabel('x[n]', 'FontSize', 12);
title(sprintf('\\kappa = %.2f, \\beta = %.2f', kappa, beta), 'FontSize', 14);
end
% Improve overall aesthetics
sgtitle('Stationary States for Different \kappa and \beta Values', 'FontSize', 16); % Super title for the figure
Mathematica plots
3 Comments
Accepted Answer
Torsten
on 18 Mar 2024
Edited: Torsten
on 18 Mar 2024
All you see in the graphs are floating point errors.
Your system of equations has solution x(i) = 0 for 1 <= i <= numBases and all three parameter constellations.
17 Comments
Torsten
on 20 Mar 2024
Edited: Torsten
on 20 Mar 2024
Because I am begginer to Matlab, when I use fsolve I replace all the following code
Yes.
and is more accurate right?
Obviously not in your case since it doesn't converge for the initialGuess you supply. Only if you set
initialGuess = 10 * ones(numBases, 1);
e.g., it converges to the solution you get with your code of Newton's method (except for the first case where at least three solutions seem to exist: 0, your solution and the solution you get with "fsolve").
So it's hard to interprete the results if you cannot apply physical reasoning to sort out solutions that make no sense.
More Answers (1)
See Also
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!