Mercurial > repos > public > sbplib

comparison +scheme/Schrodinger2dCurve.m @ 500:83734c26b8e3 feature/quantumTriangles

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
Small changes
author Ylva Rydin <ylva.rydin@telia.com>
date Thu, 20 Apr 2017 07:55:09 +0200
parents f1465e6aeb26
children ba92b27da5a0
comparison
equal deleted inserted replaced
499:f1465e6aeb26 500:83734c26b8e3
147 x_u = obj.Du*x; 147 x_u = obj.Du*x;
148 x_v = obj.Dv*x; 148 x_v = obj.Dv*x;
149 y_u = obj.Du*y; 149 y_u = obj.Du*y;
150 y_v = obj.Dv*y; 150 y_v = obj.Dv*y;
151 151
152 J = x_u.*y_v - x_v.*y_u; 152 J = (x_u.*y_v - x_v.*y_u);
153 a11 = 1./J.* (x_v.^2 + y_v.^2); 153
154 a12 = -1./J .* (x_u.*x_v + y_u.*y_v); 154 Ji = spdiags(1./J, 0, obj.m_tot, obj.m_tot);
155 a22 = 1./J .* (x_u.^2 + y_u.^2); 155 obj.Ji = Ji;
156
157 a11 = Ji* (x_v.^2 + y_v.^2);
158 a12 = -Ji* (x_u.*x_v + y_u.*y_v);
159 a22 = Ji* (x_u.^2 + y_u.^2);
156 160
157 obj.a11 = a11; 161 obj.a11 = a11;
158 obj.a12 = a12; 162 obj.a12 = a12;
159 obj.a22 = a22; 163 obj.a22 = a22;
160 164
174 D = obj.D2_v(a22(obj.ind(i,:))); 178 D = obj.D2_v(a22(obj.ind(i,:)));
175 p = obj.ind(i,:); 179 p = obj.ind(i,:);
176 obj.DVV(p,p) = D; 180 obj.DVV(p,p) = D;
177 end 181 end
178 182
179 Ji = spdiags(1./J, 0, obj.m_tot, obj.m_tot);
180 obj.Ji = Ji;
181 obj.g_1 = x_v.*y_tau-x_tau.*y_v; 183 obj.g_1 = x_v.*y_tau-x_tau.*y_v;
182 obj.g_2 = x_tau.*y_u - y_tau.*x_u; 184 obj.g_2 = x_tau.*y_u - y_tau.*x_u;
183 185
184 obj.b1 = spdiags(obj.g_1, 0, obj.m_tot, obj.m_tot); 186 obj.b1 = spdiags(obj.g_1, 0, obj.m_tot, obj.m_tot);
185 obj.b2 = spdiags(obj.g_2, 0, obj.m_tot, obj.m_tot); 187 obj.b2 = spdiags(obj.g_2, 0, obj.m_tot, obj.m_tot);
210 tau2 = @(t) (-1*s*a(t) - abs(a(t)))/4; 212 tau2 = @(t) (-1*s*a(t) - abs(a(t)))/4;
211 213
212 penalty_parameter_1 = @(t) 1*1i*halfnorm_inv_n*halfnorm_inv_t*F(t)*e'*halfnorm_t*e; 214 penalty_parameter_1 = @(t) 1*1i*halfnorm_inv_n*halfnorm_inv_t*F(t)*e'*halfnorm_t*e;
213 penalty_parameter_2 = @(t) halfnorm_inv_n*e*tau2(t); 215 penalty_parameter_2 = @(t) halfnorm_inv_n*e*tau2(t);
214 216
215 closure = @(t) obj.Ji*obj.c^2 * penalty_parameter_1(t)*e' + obj.Ji* penalty_parameter_2(t)*e'; 217 closure = @(t) (obj.Ji*obj.c^2 * penalty_parameter_1(t)*e' + obj.Ji* penalty_parameter_2(t)*e');
216 penalty = @(t) -obj.Ji*obj.c^2 * penalty_parameter_1(t)*e'- obj.Ji*penalty_parameter_2(t)*e'; 218 penalty = @(t) -obj.Ji*obj.c^2 * penalty_parameter_1(t)*e'- obj.Ji*penalty_parameter_2(t)*e';
217 219
218 end 220 end
219 221
220 function [closure, penalty] = interface(obj,boundary,neighbour_scheme,neighbour_boundary) 222 function [closure, penalty] = interface(obj,boundary,neighbour_scheme,neighbour_boundary)