From c73a80474dfe71a6fe38f779d0d761fbef535774 Mon Sep 17 00:00:00 2001
From: Anthony Royer <anthony.royer@uliege.be>
Date: Mon, 26 Sep 2022 18:09:30 +0200
Subject: [PATCH] Working PML with hetero

---
 examples/helmholtz/crossPoints/Subproblem2D.cpp |  2 +-
 examples/navier/crossPoints/Subproblem2D.cpp    | 12 ++++++------
 examples/navier/crossPoints/ddm2D.cpp           |  6 +++---
 3 files changed, 10 insertions(+), 10 deletions(-)

diff --git a/examples/helmholtz/crossPoints/Subproblem2D.cpp b/examples/helmholtz/crossPoints/Subproblem2D.cpp
index d481e858..1b315dca 100755
--- a/examples/helmholtz/crossPoints/Subproblem2D.cpp
+++ b/examples/helmholtz/crossPoints/Subproblem2D.cpp
@@ -389,7 +389,7 @@ namespace D2 {
       sigma = sigmaStar * distSigma * distSigma / (_size * _size);
     }
     
-    gmshfem::function::ScalarFunction< std::complex< double > > kMod = 1.;// + im * sigma;
+    gmshfem::function::ScalarFunction< std::complex< double > > kMod = 1. + im * sigma / kappaMod;
     if(orientation() == "E" || orientation() == "W") {
       D = gmshfem::function::tensorDiag< std::complex< double > >(1./kMod, kMod, kMod);
       E = kMod;
diff --git a/examples/navier/crossPoints/Subproblem2D.cpp b/examples/navier/crossPoints/Subproblem2D.cpp
index 3bb4864b..4e9f06c3 100755
--- a/examples/navier/crossPoints/Subproblem2D.cpp
+++ b/examples/navier/crossPoints/Subproblem2D.cpp
@@ -314,7 +314,7 @@ namespace D2 {
       sigma = sigmaStar * distSigma * distSigma / (_size * _size);
     }
     
-    gmshfem::function::ScalarFunction< std::complex< double > > gamma = 1. + im * sigma;
+    gmshfem::function::ScalarFunction< std::complex< double > > gamma = 1. + im * sigma / kappaP;
     gmshfem::function::ScalarFunction< std::complex< double > > invGamma = 1./gamma;
     if(orientation() == "E" || orientation() == "W") {
       gmshfem::function::TensorFunction< std::complex< double > > C_1 = gmshfem::function::tensor< std::complex< double > >(invGamma, invGamma, invGamma,
@@ -385,7 +385,7 @@ namespace D2 {
       sigma = sigmaStar * distSigma * distSigma / (_size * _size);
     }
     
-    gmshfem::function::ScalarFunction< std::complex< double > > kMod = 1. + im * sigma;
+    gmshfem::function::ScalarFunction< std::complex< double > > kMod = 1. + im * sigma / kappaP;
     if(orientation() == "E" || orientation() == "W") {
       D = gmshfem::function::tensorDiag< std::complex< double > >(1./kMod, kMod, kMod);
       E = kMod;
@@ -467,8 +467,8 @@ namespace D2 {
         for(int k = 0; k < 3; ++k) {
           for(int l = 0; l < 3; ++l) {
             if(i == k && j == l) {
-              tmpX(i,j)(k,l) += 1.;
-              tmpY(i,j)(k,l) += 1.;
+              tmpX(i,j)(k,l) = 1.;
+              tmpY(i,j)(k,l) = 1.;
             }
             else {
               tmpX(i,j)(k,l) = 0.;
@@ -914,8 +914,8 @@ namespace D2 {
         for(int k = 0; k < 3; ++k) {
           for(int l = 0; l < 3; ++l) {
             if(i == k && j == l) {
-              tmpX(i,j)(k,l) += 1.;
-              tmpY(i,j)(k,l) += 1.;
+              tmpX(i,j)(k,l) = 1.;
+              tmpY(i,j)(k,l) = 1.;
             }
             else {
               tmpX(i,j)(k,l) = 0.;
diff --git a/examples/navier/crossPoints/ddm2D.cpp b/examples/navier/crossPoints/ddm2D.cpp
index 134944e1..5839f1a5 100755
--- a/examples/navier/crossPoints/ddm2D.cpp
+++ b/examples/navier/crossPoints/ddm2D.cpp
@@ -440,7 +440,7 @@ namespace D2 {
         }
       }
       
-      gmshfem::function::TensorFunction< std::complex< double > > Cxx = gmshfem::function::tensor< std::complex< double > >(p1, 0., 0.,  0., p3, 0.,  0., 0., p3);
+      gmshfem::function::TensorFunction< std::complex< double > > Cxx = gmshfem::function::tensor< std::complex< double > >(p1, 0., 0.,  0., p3, 0.,  0., 0., 0.);
       gmshfem::function::TensorFunction< std::complex< double > > Cxy = gmshfem::function::tensor< std::complex< double > >(0., p2, 0.,  p2, 0., 0.,  0., 0., 0.);
       gmshfem::function::TensorFunction< std::complex< double > > Cyx = gmshfem::function::tensor< std::complex< double > >(0., p2, 0.,  p2, 0., 0.,  0., 0., 0.);
       gmshfem::function::TensorFunction< std::complex< double > > Cyy = gmshfem::function::tensor< std::complex< double > >(p3, 0., 0.,  0., p1, 0.,  0., 0., 0.);
@@ -892,8 +892,8 @@ namespace D2 {
         unsigned int index = i * nDomY + j;
         if(gmshddm::mpi::isItMySubdomain(index)) {
           if(saveU) {
-            gmshfem::post::save(xComp(u(index)), omega(index), "ux_" + std::to_string(index), "msh", "/scratch/ulg/ace/aroyer/pml/");
-            gmshfem::post::save(yComp(u(index)), omega(index), "uy_" + std::to_string(index), "msh", "/scratch/ulg/ace/aroyer/pml/");
+            gmshfem::post::save(xComp(u(index)), omega(index), "ux_" + std::to_string(index), "pos", "/scratch/ulg/ace/aroyer/pml/");
+            gmshfem::post::save(yComp(u(index)), omega(index), "uy_" + std::to_string(index), "pos", "/scratch/ulg/ace/aroyer/pml/");
             
             //gmshfem::post::save(tr(grad(u(index))), omega(index), "uP_" + std::to_string(index), "msh");
           }
-- 
GitLab