Mercurial > repos > public > sbplib
comparison +scheme/Schrodinger2dCurve.m @ 515:ba92b27da5a0 feature/quantumTriangles
Tried to make skewsym.
Values left on the diagonal when time-dep!
author | Ylva Rydin <ylva.rydin@telia.com> |
---|---|
date | Fri, 07 Jul 2017 14:22:02 +0200 |
parents | 83734c26b8e3 |
children | afff85574ddb |
comparison
equal
deleted
inserted
replaced
514:32a24485f3e8 | 515:ba92b27da5a0 |
---|---|
118 obj.variable_update(0); | 118 obj.variable_update(0); |
119 end | 119 end |
120 | 120 |
121 function [D] = d_fun(obj,t) | 121 function [D] = d_fun(obj,t) |
122 % obj.update_vairables(t); In driscretization? | 122 % obj.update_vairables(t); In driscretization? |
123 D = obj.Ji*(-1/2*(obj.b1*obj.Du-obj.b1_u+obj.Du*obj.b1) - 1/2*(obj.b2*obj.Dv - obj.b2_v +obj.Dv*obj.b2) + 1i*obj.c^2*(obj.DUU + obj.DUV + obj.DVU + obj.DVV)); | 123 % D = obj.Ji*(-1/2*(obj.b1*obj.Du-obj.b1_u+obj.Du*obj.b1) - |
124 | 124 % 1/2*(obj.b2*obj.Dv - obj.b2_v +obj.Dv*obj.b2) + |
125 % 1i*obj.c^2*(obj.DUU + obj.DUV + obj.DVU + obj.DVV)); (ols | |
126 % not skew sym disc | |
127 | |
128 D = sqrt(obj.Ji)*(-1/2*(obj.b1*obj.Du + obj.Du*obj.b1) - 1/2*(obj.b2*obj.Dv + obj.Dv*obj.b2) + 1i*obj.c^2*(obj.DUU + obj.DUV + obj.DVU + obj.DVV))*sqrt(obj.Ji); | |
125 end | 129 end |
126 | 130 |
127 | 131 |
128 function [D ]= variable_update(obj,t) | 132 function [D ]= variable_update(obj,t) |
129 % Metric derivatives | 133 % Metric derivatives |
212 tau2 = @(t) (-1*s*a(t) - abs(a(t)))/4; | 216 tau2 = @(t) (-1*s*a(t) - abs(a(t)))/4; |
213 | 217 |
214 penalty_parameter_1 = @(t) 1*1i*halfnorm_inv_n*halfnorm_inv_t*F(t)*e'*halfnorm_t*e; | 218 penalty_parameter_1 = @(t) 1*1i*halfnorm_inv_n*halfnorm_inv_t*F(t)*e'*halfnorm_t*e; |
215 penalty_parameter_2 = @(t) halfnorm_inv_n*e*tau2(t); | 219 penalty_parameter_2 = @(t) halfnorm_inv_n*e*tau2(t); |
216 | 220 |
217 closure = @(t) (obj.Ji*obj.c^2 * penalty_parameter_1(t)*e' + obj.Ji* penalty_parameter_2(t)*e'); | 221 closure = @(t) sqrt(obj.Ji)*(obj.c^2 * penalty_parameter_1(t)*e' + penalty_parameter_2(t)*e')*sqrt(obj.Ji); |
218 penalty = @(t) -obj.Ji*obj.c^2 * penalty_parameter_1(t)*e'- obj.Ji*penalty_parameter_2(t)*e'; | 222 penalty = @(t) -sqrt(obj.Ji)*(obj.c^2 * penalty_parameter_1(t)*e' - penalty_parameter_2(t)*e')*sqrt(obj.Ji); |
219 | 223 |
220 end | 224 end |
221 | 225 |
222 function [closure, penalty] = interface(obj,boundary,neighbour_scheme,neighbour_boundary) | 226 function [closure, penalty] = interface(obj,boundary,neighbour_scheme,neighbour_boundary) |
223 % u denotes the solution in the own domain | 227 % u denotes the solution in the own domain |
243 penalty_parameter_1 = @(t) halfnorm_inv_u_n*(tau*halfnorm_inv_u_t*F_u(t)*e_u'*halfnorm_u_t*e_u); | 247 penalty_parameter_1 = @(t) halfnorm_inv_u_n*(tau*halfnorm_inv_u_t*F_u(t)*e_u'*halfnorm_u_t*e_u); |
244 penalty_parameter_2 = @(t) halfnorm_inv_u_n * e_u * (sig ); | 248 penalty_parameter_2 = @(t) halfnorm_inv_u_n * e_u * (sig ); |
245 penalty_parameter_3 = @(t) halfnorm_inv_u_n * e_u * (gamm(t) ); | 249 penalty_parameter_3 = @(t) halfnorm_inv_u_n * e_u * (gamm(t) ); |
246 | 250 |
247 | 251 |
248 closure =@(t) obj.Ji*obj.c^2 * ( penalty_parameter_1(t)*e_u' + penalty_parameter_2(t)*F_u(t)' + penalty_parameter_3(t)*e_u'); | 252 closure =@(t) sqrt(obj.Ji)*obj.c^2 * ( penalty_parameter_1(t)*e_u' + penalty_parameter_2(t)*F_u(t)' + penalty_parameter_3(t)*e_u')*sqrt(obj.Ji); |
249 penalty =@(t) obj.Ji*obj.c^2 * ( -penalty_parameter_1(t)*e_v' - penalty_parameter_2(t)*F_v(t)' - penalty_parameter_3(t)*e_v'); | 253 penalty =@(t) sqrt(obj.Ji)*obj.c^2 * ( -penalty_parameter_1(t)*e_v' - penalty_parameter_2(t)*F_v(t)' - penalty_parameter_3(t)*e_v')*sqrt(obj.Ji); |
250 end | 254 end |
251 | 255 |
252 | 256 |
253 function [e, d_n, d_t, coeff_t,coeff_n, s, halfnorm_inv_n, halfnorm_inv_t, halfnorm_t,g] = get_boundary_ops(obj, boundary) | 257 function [e, d_n, d_t, coeff_t,coeff_n, s, halfnorm_inv_n, halfnorm_inv_t, halfnorm_t,g] = get_boundary_ops(obj, boundary) |
254 | 258 |