Trying to do central difference method

22 vues (au cours des 30 derniers jours)
Joe Bennett
Joe Bennett le 13 Nov 2020
Commenté : Joe Bennett le 13 Nov 2020
%QUESTION3 Implements the central differnce method to find the gradient of
%y=e^x-2sin^2(2x) at the point where x=1.
clear all
clc
format long
x=1; % defining x
y=1; % defining y, (y=x)
h=0.1 % first value of h
x=x+h; % getting x+h within f(x)
y=y-h; % getting x-h within f(x)
f = @(x) exp(x)-2*(sin(2*x)).^2; % defining f(x)
g = @(y) exp(y)-2*(sin(2*y)).^2; % defining f(y) to later use as f(x-h)
Grad=(f(x)-g(y))/2*h % Central difference method, Grad = Gradient
for n=1:2 % 2 iterations after h=0.1
h=h/10 % new value for h
x=x+h; % to replace f(x) with f(x+h)
y=y-h; % to replace f(x) with f(x-h)
Grad=(f(x)-g(y))/2*h % Central difference method
end
Here is my code, but my output is incorrect. I did the calculations on paper and they're quite different to my answers. My first value of "Grad" is correct but off by a factor of 1/100. My other two answers are completely wrong.

Réponse acceptée

John D'Errico
John D'Errico le 13 Nov 2020
Modifié(e) : John D'Errico le 13 Nov 2020
You need to understand this:
6/2*3
ans = 9
Why does it return 9, and not 1?
I told matlab to divide 6 by 2, and then to multiply the result by 3. So the answer is 9. Do you see why the following result is different?
6/(2*3)
ans = 1
Why am I pointing this out? What did you write?
Grad=(f(x)-g(y))/2*h % Central differnce method, Grad = Gradient
h needs to be in the DENOMINATOR. You want to divde by h, yet you multiplied by h. Ergo, the factor of 100, because h was 0.1.
  1 commentaire
Joe Bennett
Joe Bennett le 13 Nov 2020
Thank you, its all sorted now.

Connectez-vous pour commenter.

Plus de réponses (1)

Alan Stevens
Alan Stevens le 13 Nov 2020
You should probably use the following to calculate the gradient
Grad = (f(1+h) - f(1-h))/(2*h);
though a much smaller value for h will be needed for an accurate result.
  1 commentaire
Joe Bennett
Joe Bennett le 13 Nov 2020
Thank you, its all sorted now.

Connectez-vous pour commenter.

Catégories

En savoir plus sur Physics dans Help Center et File Exchange

Community Treasure Hunt

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

Start Hunting!

Translated by