Commit 91fb8e67 authored by Dany Dumont's avatar Dany Dumont
Browse files

Ajout de RK.m

parent 02754096
% RK - Transporte des particules avec RK4.
%
% Syntax:
%
% Inputs:
%
% Outputs:
%
% Example:
% Line 1 of example
% Line 2 of example
% Line 3 of example
%
% Other m-files required: loadheadernc.m, loadromnc.m
% Subfunctions: none
% MAT-files required: none
%
% See also: OTHER_FUNCTION_NAME1, OTHER_FUNCTION_NAME2
% Detailed description:
% Projet dans le cadre du contrat MPO-IMSER 2013-2014.
%
% Author: Simon St-Onge Drouin
% ISMER-UQAR
% email: stonsi01@uqar.ca
% Website: http://www.
% Janvier 2014; Last revision:
%------------- BEGIN CODE --------------opengl neverselect
function [k,iii] = RK(u1,u2,v1,v2,umid,vmid,k,iii,h,dx,kcomp,icomp)
% $$$ k_tmp = ( 1:length(k(end,:)) ) * NaN;
% $$$ iii_tmp = ( 1:length(k(end,:)) ) * NaN;
k_tmp = [];
iii_tmp = [];
% Vectorisé
FU1 = TriScatteredInterp(kcomp',icomp',u1(:));
clear u1
FV1 = TriScatteredInterp(kcomp',icomp',v1(:));
clear v1
FUmid = TriScatteredInterp(kcomp',icomp',umid(:));
clear umid
FVmid = TriScatteredInterp(kcomp',icomp',vmid(:));
clear vmid
FU2 = TriScatteredInterp(kcomp',icomp',u2(:));
clear u2
FV2 = TriScatteredInterp(kcomp',icomp',v2(:));
clear v2
% Calculer les 4 coefficients de Runge-Kutta
U1 = FU1(k(end,:),iii(end,:)); % K1
V1 = FV1(k(end,:),iii(end,:));
U1(isnan(U1)) = 0;
V1(isnan(V1)) = 0;
k1x = h/dx*[U1];
k1y = h/dx*[V1];
Umid = FUmid(k(end,:)+(k1x)/2,iii(end,:)+(k1y)/2); % K2
Vmid = FVmid(k(end,:)+(k1x)/2,iii(end,:)+(k1y)/2);
Umid(isnan(Umid)) = 0;
Vmid(isnan(Vmid)) = 0;
k2x = h/dx*[Umid];
k2y = h/dx*[Vmid];
Umid = FUmid(k(end,:)+(k2x)/2,iii(end,:)+(k2y)/2); % K3
Vmid = FVmid(k(end,:)+(k2x)/2,iii(end,:)+(k2y)/2);
Umid(isnan(Umid)) = 0;
Vmid(isnan(Vmid)) = 0;
k3x = h/dx*[Umid];
k3y = h/dx*[Vmid];
U2 = FU2(k(end,:)+(k3x),iii(end,:)+(k3y)); % K4
V2 = FV2(k(end,:)+(k3x),iii(end,:)+(k3y));
U2(isnan(U2)) = 0;
V2(isnan(V2)) = 0;
k4x = h/dx*[U2];
k4y = h/dx*[V2];
% Transporter les particules
kk_tmp = k(end,:) + (k1x + 2*k2x + 2*k3x + k4x)/6;
ii_tmp = iii(end,:) + (k1y + 2*k2y + 2*k3y + k4y)/6;
k_tmp = [ k_tmp; kk_tmp ];
iii_tmp = [ iii_tmp; ii_tmp ];
k = [k;k_tmp];
iii = [iii;iii_tmp];
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment