changeset 1252:8fc2f9a4c882 feature/poroelastic

Add (temporary) fix for e-s , w-n, and x-x couplings
author Martin Almquist <malmquist@stanford.edu>
date Tue, 07 Jan 2020 12:49:09 -0800
parents 02dfe3a56742
children 89dad61cad22
files +scheme/Elastic2dVariableAnisotropic.m
diffstat 1 files changed, 32 insertions(+), 1 deletions(-) [+]
line wrap: on
line diff
--- a/+scheme/Elastic2dVariableAnisotropic.m	Sat Nov 16 14:26:06 2019 -0800
+++ b/+scheme/Elastic2dVariableAnisotropic.m	Tue Jan 07 12:49:09 2020 -0800
@@ -505,6 +505,31 @@
             C_v = v.C;
             m_tot_v = v.grid.N();
 
+            % Fix {'e', 's'}, {'w', 'n'}, and {'x','x'} couplings
+            flipFlag = false;
+            e_v_flip = e_v;
+            if (strcmp(boundary,'s') && strcmp(neighbour_boundary,'e')) || ...
+               (strcmp(boundary,'e') && strcmp(neighbour_boundary,'s')) || ...
+               (strcmp(boundary,'w') && strcmp(neighbour_boundary,'n')) || ...
+               (strcmp(boundary,'n') && strcmp(neighbour_boundary,'w')) || ...
+               (strcmp(boundary,'s') && strcmp(neighbour_boundary,'s')) || ...
+               (strcmp(boundary,'n') && strcmp(neighbour_boundary,'n')) || ...
+               (strcmp(boundary,'w') && strcmp(neighbour_boundary,'w')) || ...
+               (strcmp(boundary,'e') && strcmp(neighbour_boundary,'e'))
+
+                flipFlag = true;
+                e_v_flip = fliplr(e_v);
+
+                t1 = tau_v(:,1:2:end-1);
+                t2 = tau_v(:,2:2:end);
+
+                t1 = fliplr(t1);
+                t2 = fliplr(t2);
+
+                tau_v(:,1:2:end-1) = t1;
+                tau_v(:,2:2:end) = t2;
+            end
+
             % Operators that are only required for own domain
             Hi      = u.Hi_kron;
             RHOi    = u.RHOi_kron;
@@ -541,10 +566,16 @@
                             Z_v = Z_v + e_v'*nu_v(l)*C_v{l,i,k,j}*nu_v(k)*e_v;
                         end
                     end
+
+                    if flipFlag
+                        Z_v = rot90(Z_v,2);
+                    end
+
                     Z = -tuning*dim*( 1/(4*h11_u)*Z_u + 1/(4*h11_v)*Z_v );
                     X = Y + Z*e_u';
                     closure = closure + E_u{i}*X'*H_gamma*e_u'*E_u{j}';
-                    penalty = penalty - E_u{i}*X'*H_gamma*e_v'*E_v{j}';
+                    penalty = penalty - E_u{i}*X'*H_gamma*e_v_flip'*E_v{j}';
+
                 end
             end