From 28d71fcf5103946adf2050094284b27b4f9d08fb Mon Sep 17 00:00:00 2001
From: Nicolas Marsic <marsic@temf.tu-darmstadt.de>
Date: Wed, 28 Apr 2021 17:46:26 +0200
Subject: [PATCH] Starting RCot TC (looks good for 2 equal subdomains)
---
GetDDM/Helmholtz.pro | 40 +++++++++++-
GetDDM/cavity2d.pro | 133 ++++++++++++++++++++++++++++++++++++++++
GetDDM/cavity2d_src.pro | 115 ++++++++++++++++++++++++++++++++++
3 files changed, 285 insertions(+), 3 deletions(-)
create mode 100644 GetDDM/cavity2d.pro
create mode 100644 GetDDM/cavity2d_src.pro
diff --git a/GetDDM/Helmholtz.pro b/GetDDM/Helmholtz.pro
index 5f33571..cda0145 100644
--- a/GetDDM/Helmholtz.pro
+++ b/GetDDM/Helmholtz.pro
@@ -124,7 +124,7 @@ FunctionSpace {
}
}
}
- If (TC_TYPE == 2)
+ If (TC_TYPE == 2 || TC_TYPE == 4)
For h In {1:NP_OSRC}
{ Name Hgrad_phi~{h}~{i}~{j}; Type Form0;
BasisFunction {
@@ -148,7 +148,7 @@ Formulation {
{ Name u~{i}; Type Local; NameOfSpace Hgrad_u~{i}; }
For jj In {0: #myD~{i}()-1}
j = myD~{i}(jj);
- If(TC_TYPE == 2)
+ If(TC_TYPE == 2 || TC_TYPE == 4)
For h In{1:NP_OSRC}
{ Name phi~{h}~{i}~{j}; Type Local;
NameOfSpace Hgrad_phi~{h}~{i}~{j}; }
@@ -226,6 +226,28 @@ Formulation {
EndFor
EndIf
+ If(TC_TYPE == 4) // RCot
+ For jj In {0: #myD~{i}()-1}
+ j = myD~{i}(jj);
+ Galerkin { [ 1/Ld~{j} * Dof{u~{i}}, {u~{i}} ];
+ In Sigma~{i}~{j}; Jacobian JSur; Integration I1; }
+
+ For h In{1:NP_OSRC}
+ Galerkin { [ +2/Ld~{j} * Dof{phi~{h}~{i}~{j}}, {u~{i}} ];
+ In Sigma~{i}~{j}; Jacobian JSur; Integration I1; }
+ Galerkin { [ -2/Ld~{j} * (1/k[])^2 * Dof{d phi~{h}~{i}~{j}}, {d u~{i}} ];
+ In Sigma~{i}~{j}; Jacobian JSur; Integration I1; }
+
+ Galerkin { [ -1 * Dof{u~{i}}, {phi~{h}~{i}~{j}} ];
+ In Sigma~{i}~{j}; Jacobian JSur; Integration I1; }
+ Galerkin { [ (1-((h*Pi)/(k[]*Ld~{j}))^2) * Dof{phi~{h}~{i}~{j}}, {phi~{h}~{i}~{j}} ];
+ In Sigma~{i}~{j}; Jacobian JSur; Integration I1; }
+ Galerkin { [ -(1/k[])^2 * Dof{d phi~{h}~{i}~{j}}, {d phi~{h}~{i}~{j}} ];
+ In Sigma~{i}~{j}; Jacobian JSur; Integration I1; }
+ EndFor
+ EndFor
+ EndIf
+
// Bayliss-Turkel absorbing boundary condition
Galerkin { [ - I[] * kInf[] * Dof{u~{i}} , {u~{i}} ];
In GammaInf~{i}; Jacobian JSur; Integration I1; }
@@ -246,7 +268,7 @@ Formulation {
{ Name u~{i}; Type Local; NameOfSpace Hgrad_u~{i}; }
{ Name g_out~{i}~{j}; Type Local;
NameOfSpace Hgrad_g_out~{i}~{j}; }
- If(TC_TYPE == 2)
+ If(TC_TYPE == 2 || TC_TYPE == 4)
For h In{1:NP_OSRC}
{ Name phi~{h}~{i}~{j}; Type Local;
NameOfSpace Hgrad_phi~{h}~{i}~{j}; }
@@ -290,6 +312,18 @@ Formulation {
Galerkin { [ 2 * I[] * kInf[] * {u~{i}}, {g_out~{i}~{j}}];
In TrBndPmlSigma~{i}~{j}; Jacobian JSur; Integration I1;}
EndIf
+
+ If(TC_TYPE == 4) // RCot
+ Galerkin { [ -2 * 1/Ld~{j} * {u~{i}}, {g_out~{i}~{j}} ];
+ In Sigma~{i}~{j}; Jacobian JSur; Integration I1; }
+ For h In{1:NP_OSRC}
+ Galerkin { [ -2 * 2/Ld~{j} * {u~{i}}, {g_out~{i}~{j}} ];
+ In Sigma~{i}~{j}; Jacobian JSur; Integration I1; }
+ Galerkin { [ -2 * 2/Ld~{j} * ((h*Pi)/(k[]*Ld~{j}))^2 * {phi~{h}~{i}~{j}}, {g_out~{i}~{j}} ];
+ In Sigma~{i}~{j}; Jacobian JSur; Integration I1; }
+ EndFor
+ EndIf
+
}
}
diff --git a/GetDDM/cavity2d.pro b/GetDDM/cavity2d.pro
new file mode 100644
index 0000000..3509f3e
--- /dev/null
+++ b/GetDDM/cavity2d.pro
@@ -0,0 +1,133 @@
+Include "cavity2d_data.geo";
+
+DefineConstant[
+ TC_TYPE = {2, Name "Input/01Transmission condition",
+ Choices {0="Order 0", 1="Order 2", 2="Pade (OSRC)", 4="RCot"}},
+ NP_OSRC = {4, Name "Input/02Number of auxiliary variables (OSRC or RCot)"},
+ PRECONDITIONER = 0,
+ IS_OPEN = {0, Name "Input/03Is cavity open?", Choices{0, 1}},
+ MODE = {-1, Name "Input/04Mode"},
+ ListOfCuts() = { {0, N_DOM-1} }
+];
+
+Include "cavity2d_src.pro";
+
+Function{
+ // Imaginary value
+ I[] = Complex[0, 1];
+
+ // Domain length
+ Ld~{0} = LX * t2Dom;
+ Ld~{1} = LX * (1-t2Dom);
+
+ // Wavenumber
+ k = WAVENUMBER;
+ k[] = k;
+
+ // ABC
+ kInf[] = k;
+ alphaBT[] = 0;
+ betaBT[] = 0;
+
+ // 0th order TC
+ kIBC[] = k[];
+
+ // 2nd order TC (Gander 2002, pp. 46-47)
+ xsimin = 0;
+ xsimax = Pi / LC;
+ deltak[] = Pi;
+ alphastar[] = I[] * ((k[]^2 - xsimin^2) * (k[]^2 - (k[]-deltak[])^2))^(1/4);
+ betastar[] = ((xsimax^2 - k[]^2) * ((k[]+deltak[])^2 - k[]^2))^(1/4);
+ a[] = - (alphastar[] * betastar[] - k[]^2) / (alphastar[] + betastar[]);
+ b[] = - 1 / (alphastar[] + betastar[]);
+
+ // Pade TC
+ kepsI = 0;
+ keps[] = k[]*(1+kepsI*I[]);
+ theta_branch = Pi/4;
+
+ // No PMLs
+ D[] = 1;
+ E[] = 1;
+}
+
+Group{
+ D() = {};
+ For idom In {0:N_DOM-1}
+ left = (idom-1)%N_DOM; // left neighbour
+ right = (idom+1)%N_DOM; // right neighbour
+
+ D() += idom;
+
+ Omega~{idom} = Region[(idom + 1)];
+ GammaD~{idom} = Region[{}];
+ GammaD0~{idom} = Region[(idom + 2*N_DOM + 2)];
+ GammaInf~{idom} = Region[{}];
+ GammaN~{idom} = Region[{}];
+
+ If(idom == 0)
+ D~{idom} = {right};
+
+ Sigma~{idom}~{right} = Region[(idom + N_DOM + 2)];
+ Sigma~{idom} = Region[{Sigma~{idom}~{right}}];
+ GammaD~{idom} += Region[(idom + N_DOM + 1)];
+
+ BndSigma~{idom}~{right} = Region[{}];
+ BndSigma~{idom} = Region[{BndSigma~{idom}~{right}}];
+ BndGammaInf~{idom}~{right} = Region[{}];
+ BndGammaInf~{idom} = Region[{BndGammaInf~{idom}~{right}}];
+
+ Pml~{idom}~{right} = Region[{}];
+ PmlD0~{idom}~{right} = Region[{}];
+ PmlInf~{idom}~{right} = Region[{}];
+ EndIf
+ If(idom == N_DOM-1)
+ D~{idom} = {left};
+
+ Sigma~{idom}~{left} = Region[(idom + N_DOM + 1)];
+ Sigma~{idom} = Region[{Sigma~{idom}~{left}}];
+
+ If(IS_OPEN)
+ GammaInf~{idom} += Region[(idom + N_DOM + 2)];
+ Else
+ GammaD0~{idom} += Region[(idom + N_DOM + 2)];
+ EndIf
+ BndSigma~{idom}~{left} = Region[{}];
+ BndSigma~{idom} = Region[{BndSigma~{idom}~{left}}];
+ BndGammaInf~{idom}~{left} = Region[{}];
+ BndGammaInf~{idom} = Region[{BndGammaInf~{idom}~{left}}];
+
+ Pml~{idom}~{left} = Region[{}];
+ PmlD0~{idom}~{left} = Region[{}];
+ PmlInf~{idom}~{left} = Region[{}];
+ EndIf
+ If(idom >= 1 && idom < N_DOM-1)
+ D~{idom} = {left, right};
+
+ Sigma~{idom}~{left} = Region[(idom + N_DOM + 1)];
+ Sigma~{idom}~{right} = Region[(idom + N_DOM + 2)];
+ Sigma~{idom} = Region[{Sigma~{idom}~{left},
+ Sigma~{idom}~{right}}];
+
+ BndSigma~{idom}~{left} = Region[{}];
+ BndSigma~{idom}~{right} = Region[{}];
+ BndSigma~{idom} = Region[{BndSigma~{idom}~{left},
+ BndSigma~{idom}~{right}}];
+ BndGammaInf~{idom}~{left} = Region[{}];
+ BndGammaInf~{idom}~{right} = Region[{}];
+ BndGammaInf~{idom} = Region[{BndGammaInf~{idom}~{left},
+ BndGammaInf~{idom}~{right}}];
+
+ Pml~{idom}~{left} = Region[{}];
+ Pml~{idom}~{right} = Region[{}];
+ PmlD0~{idom}~{left} = Region[{}];
+ PmlD0~{idom}~{right} = Region[{}];
+ PmlInf~{idom}~{left} = Region[{}];
+ PmlInf~{idom}~{right} = Region[{}];
+ EndIf
+ EndFor
+}
+
+Include "Decomposition.pro";
+Include "Helmholtz.pro";
+Include "Schwarz.pro";
diff --git a/GetDDM/cavity2d_src.pro b/GetDDM/cavity2d_src.pro
new file mode 100644
index 0000000..8cb4a59
--- /dev/null
+++ b/GetDDM/cavity2d_src.pro
@@ -0,0 +1,115 @@
+Function{
+ If(MODE > 0)
+ ky = MODE * Pi / LY;
+ uinc[] = Complex[Sin[ky * Y[]], 0];
+ Else
+ uinc[] = Complex[Sin[ 1 * Pi/LY * Y[]], 0] +
+ Complex[Sin[ 2 * Pi/LY * Y[]], 0] +
+ Complex[Sin[ 3 * Pi/LY * Y[]], 0] +
+ Complex[Sin[ 4 * Pi/LY * Y[]], 0] +
+ Complex[Sin[ 5 * Pi/LY * Y[]], 0] +
+ Complex[Sin[ 6 * Pi/LY * Y[]], 0] +
+ Complex[Sin[ 7 * Pi/LY * Y[]], 0] +
+ Complex[Sin[ 8 * Pi/LY * Y[]], 0] +
+ Complex[Sin[ 9 * Pi/LY * Y[]], 0] +
+
+ Complex[Sin[10 * Pi/LY * Y[]], 0] +
+ Complex[Sin[11 * Pi/LY * Y[]], 0] +
+ Complex[Sin[12 * Pi/LY * Y[]], 0] +
+ Complex[Sin[13 * Pi/LY * Y[]], 0] +
+ Complex[Sin[14 * Pi/LY * Y[]], 0] +
+ Complex[Sin[15 * Pi/LY * Y[]], 0] +
+ Complex[Sin[16 * Pi/LY * Y[]], 0] +
+ Complex[Sin[17 * Pi/LY * Y[]], 0] +
+ Complex[Sin[18 * Pi/LY * Y[]], 0] +
+ Complex[Sin[19 * Pi/LY * Y[]], 0] +
+
+ Complex[Sin[20 * Pi/LY * Y[]], 0] +
+ Complex[Sin[21 * Pi/LY * Y[]], 0] +
+ Complex[Sin[22 * Pi/LY * Y[]], 0] +
+ Complex[Sin[23 * Pi/LY * Y[]], 0] +
+ Complex[Sin[24 * Pi/LY * Y[]], 0] +
+ Complex[Sin[25 * Pi/LY * Y[]], 0] +
+ Complex[Sin[26 * Pi/LY * Y[]], 0] +
+ Complex[Sin[27 * Pi/LY * Y[]], 0] +
+ Complex[Sin[28 * Pi/LY * Y[]], 0] +
+ Complex[Sin[29 * Pi/LY * Y[]], 0] +
+
+ Complex[Sin[30 * Pi/LY * Y[]], 0] +
+ Complex[Sin[31 * Pi/LY * Y[]], 0] +
+ Complex[Sin[32 * Pi/LY * Y[]], 0] +
+ Complex[Sin[33 * Pi/LY * Y[]], 0] +
+ Complex[Sin[34 * Pi/LY * Y[]], 0] +
+ Complex[Sin[35 * Pi/LY * Y[]], 0] +
+ Complex[Sin[36 * Pi/LY * Y[]], 0] +
+ Complex[Sin[37 * Pi/LY * Y[]], 0] +
+ Complex[Sin[38 * Pi/LY * Y[]], 0] +
+ Complex[Sin[39 * Pi/LY * Y[]], 0] +
+
+ Complex[Sin[40 * Pi/LY * Y[]], 0] +
+ Complex[Sin[41 * Pi/LY * Y[]], 0] +
+ Complex[Sin[42 * Pi/LY * Y[]], 0] +
+ Complex[Sin[43 * Pi/LY * Y[]], 0] +
+ Complex[Sin[44 * Pi/LY * Y[]], 0] +
+ Complex[Sin[45 * Pi/LY * Y[]], 0] +
+ Complex[Sin[46 * Pi/LY * Y[]], 0] +
+ Complex[Sin[47 * Pi/LY * Y[]], 0] +
+ Complex[Sin[48 * Pi/LY * Y[]], 0] +
+ Complex[Sin[49 * Pi/LY * Y[]], 0] +
+
+ Complex[Sin[50 * Pi/LY * Y[]], 0] +
+ Complex[Sin[51 * Pi/LY * Y[]], 0] +
+ Complex[Sin[52 * Pi/LY * Y[]], 0] +
+ Complex[Sin[53 * Pi/LY * Y[]], 0] +
+ Complex[Sin[54 * Pi/LY * Y[]], 0] +
+ Complex[Sin[55 * Pi/LY * Y[]], 0] +
+ Complex[Sin[56 * Pi/LY * Y[]], 0] +
+ Complex[Sin[57 * Pi/LY * Y[]], 0] +
+ Complex[Sin[58 * Pi/LY * Y[]], 0] +
+ Complex[Sin[59 * Pi/LY * Y[]], 0] +
+
+ Complex[Sin[60 * Pi/LY * Y[]], 0] +
+ Complex[Sin[61 * Pi/LY * Y[]], 0] +
+ Complex[Sin[62 * Pi/LY * Y[]], 0] +
+ Complex[Sin[63 * Pi/LY * Y[]], 0] +
+ Complex[Sin[64 * Pi/LY * Y[]], 0] +
+ Complex[Sin[65 * Pi/LY * Y[]], 0] +
+ Complex[Sin[66 * Pi/LY * Y[]], 0] +
+ Complex[Sin[67 * Pi/LY * Y[]], 0] +
+ Complex[Sin[68 * Pi/LY * Y[]], 0] +
+ Complex[Sin[69 * Pi/LY * Y[]], 0] +
+
+ Complex[Sin[70 * Pi/LY * Y[]], 0] +
+ Complex[Sin[71 * Pi/LY * Y[]], 0] +
+ Complex[Sin[72 * Pi/LY * Y[]], 0] +
+ Complex[Sin[73 * Pi/LY * Y[]], 0] +
+ Complex[Sin[74 * Pi/LY * Y[]], 0] +
+ Complex[Sin[75 * Pi/LY * Y[]], 0] +
+ Complex[Sin[76 * Pi/LY * Y[]], 0] +
+ Complex[Sin[77 * Pi/LY * Y[]], 0] +
+ Complex[Sin[78 * Pi/LY * Y[]], 0] +
+ Complex[Sin[79 * Pi/LY * Y[]], 0] +
+
+ Complex[Sin[80 * Pi/LY * Y[]], 0] +
+ Complex[Sin[81 * Pi/LY * Y[]], 0] +
+ Complex[Sin[82 * Pi/LY * Y[]], 0] +
+ Complex[Sin[83 * Pi/LY * Y[]], 0] +
+ Complex[Sin[84 * Pi/LY * Y[]], 0] +
+ Complex[Sin[85 * Pi/LY * Y[]], 0] +
+ Complex[Sin[86 * Pi/LY * Y[]], 0] +
+ Complex[Sin[87 * Pi/LY * Y[]], 0] +
+ Complex[Sin[88 * Pi/LY * Y[]], 0] +
+ Complex[Sin[89 * Pi/LY * Y[]], 0] +
+
+ Complex[Sin[90 * Pi/LY * Y[]], 0] +
+ Complex[Sin[91 * Pi/LY * Y[]], 0] +
+ Complex[Sin[92 * Pi/LY * Y[]], 0] +
+ Complex[Sin[93 * Pi/LY * Y[]], 0] +
+ Complex[Sin[94 * Pi/LY * Y[]], 0] +
+ Complex[Sin[95 * Pi/LY * Y[]], 0] +
+ Complex[Sin[96 * Pi/LY * Y[]], 0] +
+ Complex[Sin[97 * Pi/LY * Y[]], 0] +
+ Complex[Sin[98 * Pi/LY * Y[]], 0] +
+ Complex[Sin[99 * Pi/LY * Y[]], 0];
+ EndIf
+}
--
GitLab