CODEWITHSUNDEEP
matlabmathstochastic
David
David

Reputation: 77

Stochastic Differential Equations (SDE) in 2 dimensions

I am working on stochastic differential equations for the first time. I am looking to simulate and solve a stochastic differential equations in two dimensions.

The model is as follows:

dp=F(t,p)dt+G(t,p)dW(t)

where:

I wrote code as follows:

function MDL=gyro_2dim(Psi,D)
% want to solve for 2-by-1 vector:
%p=[theta;phi];
%drift function
F=@(t,theta,phi)  [sinth(theta)+Psi.*cos(phi)-D.*cot(theta);Psi.*cot(theta).*sin(phi)];
%diffusion function
G=@(t,theta,phi) [D 0;0 D./sin(theta)];
MDL=sde(F,G)
end

Then I call the function with the following script:

params.t0   = 0;               % start time of simulation
params.tend = 20;              % end time
params.dt =0.1;                % time increment
D=0.1;
nPeriods=10; % # of simulated observations
Psi=1;
MDL=gyro_2dim(Psi,D);
[S,T,Z]=simulate(MDL, nPeriods,'DeltaTime',params.dt);
plot(T,S)

When I run the code, I receive this error message:

Drift rate invalid at initial conditions or inconsistent model dimensions.

Any idea how to fix this error?

Upvotes: 3

Views: 1162

Answers (1)

horchler
horchler

Reputation: 18484

From the documentation for sde:

User-defined drift-rate function, denoted by F. DriftRate is a function that returns an NVARS-by-1 drift-rate vector when called with two inputs:
- A real-valued scalar observation time t.
- An NVARS-by-1 state vector Xt.

A similar specification is provided for the Diffusion function. However, you're passing in the elements of your state vector as scalars and thus have three, rather than two, inputs. You can try changing your model creation function to:

function MDL=gyro_2dim(Psi,D)
% State vector: p = [theta;phi];
F = @(t,p)[sin(p(1))+Psi.*cos(p(2))-D.*cot(p(1));
           Psi.*cot(p(1)).*sin(p(2))];            % Drift
G = @(t,p)[D 0;
           0 D./sin(p(1))];                       % Diffusion
MDL = sde(F,G);
MDL.StartTime = 0;   % Set initial time
MDL.StartState = ... % Set initial conditions

I also changed sinth(theta) to sin(p(1)) as there is no sinth function. I can't test this as I don't have the Financial toolbox (few do).

Upvotes: 3

Related Questions