Mercurial > repos > public > sbplib

annotate +time/+expint/Magnus_4.m @ 512:4ef2d2a493f1 feature/quantumTriangles

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
add H_xyz to hypsyst to compute H-norm
author Ylva Rydin <ylva.rydin@telia.com>
date Mon, 26 Jun 2017 19:23:19 +0200
parents 71908bbce2e1
children bc39bb984d88
Ignore whitespace changes - Everywhere: Within whitespace: At end of lines:
rev   line source
510
71908bbce2e1 tried to add 4th order magnus integrator which does not work
Ylva Rydin <ylva.rydin@telia.com>
parents:
diff changeset
1 % Takes one time step of size k using a fourth order magnus integrator
71908bbce2e1 tried to add 4th order magnus integrator which does not work
Ylva Rydin <ylva.rydin@telia.com>
parents:
diff changeset
2 % starting from v_0 and where the function F(v,t) gives the
71908bbce2e1 tried to add 4th order magnus integrator which does not work
Ylva Rydin <ylva.rydin@telia.com>
parents:
diff changeset
3 % time derivatives.
71908bbce2e1 tried to add 4th order magnus integrator which does not work
Ylva Rydin <ylva.rydin@telia.com>
parents:
diff changeset
4 function v = Magnus_4(v,D, t , k)
71908bbce2e1 tried to add 4th order magnus integrator which does not work
Ylva Rydin <ylva.rydin@telia.com>
parents:
diff changeset
5
71908bbce2e1 tried to add 4th order magnus integrator which does not work
Ylva Rydin <ylva.rydin@telia.com>
parents:
diff changeset
6 if isa(D,'function_handle')
71908bbce2e1 tried to add 4th order magnus integrator which does not work
Ylva Rydin <ylva.rydin@telia.com>
parents:
diff changeset
7 % v = krylov(k*D(t +k/2*t),v);
71908bbce2e1 tried to add 4th order magnus integrator which does not work
Ylva Rydin <ylva.rydin@telia.com>
parents:
diff changeset
8 c1 = 1/2 - sqrt(3)/6;
71908bbce2e1 tried to add 4th order magnus integrator which does not work
Ylva Rydin <ylva.rydin@telia.com>
parents:
diff changeset
9 c2 = 1/2 + sqrt(3)/6;
71908bbce2e1 tried to add 4th order magnus integrator which does not work
Ylva Rydin <ylva.rydin@telia.com>
parents:
diff changeset
10
71908bbce2e1 tried to add 4th order magnus integrator which does not work
Ylva Rydin <ylva.rydin@telia.com>
parents:
diff changeset
11 A1 = D(t +c1*k);
71908bbce2e1 tried to add 4th order magnus integrator which does not work
Ylva Rydin <ylva.rydin@telia.com>
parents:
diff changeset
12 A2 = D(t + c2*k);
71908bbce2e1 tried to add 4th order magnus integrator which does not work
Ylva Rydin <ylva.rydin@telia.com>
parents:
diff changeset
13 Omega = k/2*(A1 + A2) + sqrt(3)*k^2/12*(A1*A2-A2*A1);
512
4ef2d2a493f1 add H_xyz to hypsyst to compute H-norm
Ylva Rydin <ylva.rydin@telia.com>
parents: 510
diff changeset
14 % v = expm(Omega)*v;
4ef2d2a493f1 add H_xyz to hypsyst to compute H-norm
Ylva Rydin <ylva.rydin@telia.com>
parents: 510
diff changeset
15 toler = 10^(-8);
4ef2d2a493f1 add H_xyz to hypsyst to compute H-norm
Ylva Rydin <ylva.rydin@telia.com>
parents: 510
diff changeset
16 v = time.expint.expm_Arnoldi(-Omega,v,k,toler,100);
510
71908bbce2e1 tried to add 4th order magnus integrator which does not work
Ylva Rydin <ylva.rydin@telia.com>
parents:
diff changeset
17 else
71908bbce2e1 tried to add 4th order magnus integrator which does not work
Ylva Rydin <ylva.rydin@telia.com>
parents:
diff changeset
18 %v = krylov(k*D,v);
71908bbce2e1 tried to add 4th order magnus integrator which does not work
Ylva Rydin <ylva.rydin@telia.com>
parents:
diff changeset
19 v = expm(k*D)*v;
71908bbce2e1 tried to add 4th order magnus integrator which does not work
Ylva Rydin <ylva.rydin@telia.com>
parents:
diff changeset
20 end
71908bbce2e1 tried to add 4th order magnus integrator which does not work
Ylva Rydin <ylva.rydin@telia.com>
parents:
diff changeset
21
71908bbce2e1 tried to add 4th order magnus integrator which does not work
Ylva Rydin <ylva.rydin@telia.com>
parents:
diff changeset
22 end