function CINEKovn181115

ed7p10 % Kommentera bort för att köra!
%eds6p9

function sauer5p1p8

% Hitta feltermen och ordning f�r formeln:
% f'(x) = (4*f(x+h) - 3*f(x) - f(x-2h))/6h

% f(x+h) = f(x) + f'(x)*h + f''(x)/2*h^2 + f'''(x)/6*h^3 %
% Taylor-utveckling

% 6h*f' = 4*(f(x) + f'(x)*h + f''(x)/2*h^2 + f'''(x)/6*h^3) + 
%        -3*f(x) + 
%    -(f(x) - 2*f'(x)*h + f''(x)/2*(-2h)^2 + f'''(x)/6*(-2h)^3)

%  6h*f' = 4*(f(x) + f'(x)*h + f''(x)/2*h^2 + f'''(x)/6*h^3) + 
%        -3*f(x) + 
%    -(f(x) - 2*f'(x)*h + 2*f''(x)*h^2 - 4*f'''(x)/3*h^3)

% 6h*f' = 0*f(x) + 6h*f' + 0*f''*h^2 + 2*f'''*h^3

% felterm = 2*f'''*h^3 / 6h = 1/3*f'''*h^2


function ed7p10
% y10 = f(A), felet = abs(df/Da*Ea)

h = 1e-3;
A = 0.2; 
dfdA = abs((y10((1+h)*A) - y10(A))/(h*A)); % (f(x+hh) - f(x))/hh

felA = 5e-3/dfdA;

y10(A) - y10(A+felA)

function yut = y10(A)

% y'' + 2xy + y^2 = 0
% u = [y y'] -> u' = [y' y''] = [y' -(2*xy + y^2)] = [u(2); -(2*x*u(1) + u(1)^2)]

uprim = @(x, u) [u(2); -(2*x*u(1) + u(1)^2)];

% A = 0.2;

u0 = [0.3 A]';

odset = odeset('RelTol', 1e-5);

[x, u] = ode45(uprim, [0.3 10], u0, odset);

yut = u(end, 1);

function eds6p9

y = @(x) sqrt((1-x).*(1+2*x)).*(x.^2 + 0.1*x+0.1);

yx = @(x) x.*y(x);

tol = 2.5e-9;

Iy = 2*integral(y, -0.5, 1, 'RelTol', tol);

Iyx = integral(yx, -0.5, 1, 'RelTol', tol);

xtp = 2*Iyx/Iy

% maxfel i xtp blir ((1+tol)*Iyx)/((1-tol)*Iy) = (1+tol)/(1-tol)*xtp och
% minfel blir (1-tol)/(1+tol)*xtp, där den första är längst ifrån xtp, om 0 < tol < 1.

% (1+tol)/((1-tol) -1 > 1-(1-tol)/((1+tol)

% och vi väljer den första, sådant att:

% (1+tol)/((1-tol) -1 < 0.5e-8 ger 

% tol < 2.5e-9

% Monte-Carlo-beräkning av Iy ovan, genom att slumpa punkter i en rektangel
% som omsluter Iy och sedan räkna vilka som hamlar i och utanför,
% vars relation till arean av rektangeln bestämmer Iy.

N = 1e6; % Antal punkter. För stora N blir Iym och Iy väldigt lika!

Ra = 1.5*1.5

p = 1.5*rand(2,N); % Slumpa punkter i rektangeln som är 1.5*1.5 stor.

p(1,:) = p(1,:)-0.5; % Justera till rätt intervall i x-led
p(2,:) = p(2,:)-0.75; % och i y-led.

M = 0;
for i = 1:N
    ypp = y(p(1,i)); % kurvans y-koord.
    yp = p(2,i); % Slumpade punktens y-koord.
    if (yp > -ypp) && (yp < ypp) % Om y-koord. för punkten är inom [-y, +y]
        M = M + 1; % räkna som träff.
    end
end

Iym = M/N*Ra
Iy



