Mercurial > repos > public > sbplib

changeset 506:2047ff9640b3 feature/quantumTriangles

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
tried to make a symmetric scheme
author Ylva Rydin <ylva.rydin@telia.com>
date Mon, 29 May 2017 13:54:16 +0200
parents 49bc2929dfb1
children 1341e8580d64
files +scheme/Schrodinger1dCurve.m
diffstat 1 files changed, 9 insertions(+), 7 deletions(-) [+]
line wrap: on
line diff
--- a/+scheme/Schrodinger1dCurve.m	Thu May 18 11:44:29 2017 +0200
+++ b/+scheme/Schrodinger1dCurve.m	Mon May 29 13:54:16 2017 +0200
@@ -17,6 +17,7 @@
         a
         a_xi
         Ji
+        J
         t_up
         x
         
@@ -85,8 +86,8 @@
         
         function [D] = d_fun(obj,t)
             obj.variable_update(t); % In driscretization?
-            D= obj.Ji*(-0.5*(obj.D1*obj.a - obj.a_xi + obj.a*obj.D1) + 1i*obj.D2(diag(obj.Ji)) + 1i*obj.V_mat);
-            
+            D = sqrt(obj.Ji)*(-0.5*(obj.D1*obj.a*sqrt(obj.Ji) + obj.a*obj.D1*sqrt(obj.Ji)) + 1i*obj.D2(diag(obj.Ji))*obj.Ji + 1i*obj.V_mat);
+          %   D= obj.Ji*(-sqrt(obj.a)*obj.D1*sqrt(obj.a) + 0.5*obj.a_xi + 1i*obj.D2(diag(obj.Ji)) + 1i*obj.V_mat);
         end
         
         
@@ -102,6 +103,7 @@
                 obj.a = sparse(diag((-ddt_x_l*( x_r - x_l) - (obj.x-x_l)*(ddt_x_r-ddt_x_l))/(x_r-x_l)));            
                 
                 obj.Ji = sparse(diag(1./(x_r - x_l + 0*obj.x)));
+                obj.J = sparse(x_r -x_l);
                 obj.a_xi = sparse(diag(-1*(ddt_x_r - ddt_x_l + 0*obj.x)));
                 obj.t_up = t;
             end
@@ -116,17 +118,17 @@
             switch type
                 % Dirichlet boundary condition
                 case {'D','d','dirichlet'}
-                    tau1 = @(t) s * 1i*obj.Ji(p,p)^2*d;
-                    tau2 = @(t) obj.Ji*(-1*s*obj.a(p,p) - abs(obj.a(p,p)))/4*e;
-                    closure = @(t) obj.Hi*tau1(t)*e' + obj.Hi*tau2(obj.a)*e';
+                    tau1 = @(t) s * 1i*obj.Ji(p,p)*d;
+                    tau2 = @(t) (1*s*obj.a(p,p))/2*e;
+                    closure = @(t)(obj.Ji*obj.Hi*tau1(t)*e' + obj.Hi*tau2(obj.a)*e');
                     
                     switch class(data)
                         case 'double'
-                            penalty = @(t) -obj.Ji*(obj.Hi*tau1*data+obj.Hi*tau2(obj.a)*data);
+                            penalty = @(t) -(obj.Hi*tau1*data+obj.Hi*tau2(obj.a)*data);
                             %                      case 'function_handle'
                             %                           penalty = @(t)-obj.Hi*tau*data(t);
                         otherwise
-                            error('Wierd data argument!')
+                            error('Weird data argument!')
                     end
                     
                     % Unknown, boundary condition