From dcfd65cc33f5fbef3a437f0e745ff889622252ac Mon Sep 17 00:00:00 2001
From: Axel Modave <axel.modave@ensta-paristech.fr>
Date: Mon, 29 Oct 2018 21:16:06 +0100
Subject: [PATCH] AcademicWaves: update + more onelab (EM not checked)
---
AcademicWaves/formulations_scalarWaves.pro | 314 ++++++++++-----------
AcademicWaves/formulations_vectorWaves.pro | 74 ++---
AcademicWaves/gaussInCuboid.dat | 17 +-
AcademicWaves/gaussInCuboid.geo | 4 +-
AcademicWaves/gaussInCuboid.pro | 48 ++--
AcademicWaves/gaussInSquare.dat | 19 +-
AcademicWaves/gaussInSquare.geo | 12 +-
AcademicWaves/gaussInSquare.pro | 51 ++--
AcademicWaves/loudspeaker.dat | 33 +--
AcademicWaves/loudspeaker.geo | 40 +--
AcademicWaves/loudspeaker.pro | 85 ++----
AcademicWaves/scattKite.dat | 20 +-
AcademicWaves/scattKite.geo | 142 +++++-----
AcademicWaves/scattKite.pro | 95 +++----
AcademicWaves/scattPlates.dat | 16 +-
AcademicWaves/scattPlates.geo | 57 ++--
AcademicWaves/scattPlates.pro | 92 +++---
AcademicWaves/scattSphere.dat | 10 +-
AcademicWaves/scattSphere.geo | 3 +-
AcademicWaves/scattSphere.pro | 8 +-
20 files changed, 476 insertions(+), 664 deletions(-)
diff --git a/AcademicWaves/formulations_scalarWaves.pro b/AcademicWaves/formulations_scalarWaves.pro
index bd052f3..c41b08f 100644
--- a/AcademicWaves/formulations_scalarWaves.pro
+++ b/AcademicWaves/formulations_scalarWaves.pro
@@ -9,258 +9,254 @@
//==============================================================================
Function {
- DefineFunction[ fNeumann, alphaBT, betaBT, pmlScal, pmlTens ] ;
- DefineVariable[ DIM, ORDER ] ;
- DefineVariable[ dt, beta, gamma, theta, tf ] ;
+ DefineFunction[ fNeumann, alphaBT, betaBT, pmlScal, pmlTens ];
+ DefineVariable[ dt, beta, gamma, theta, tf, FREQ, k ];
}
Group{
- DefineGroup[ Omega, OmegaPml, GammaD, GammaN, GammaR, GammaBT ] ;
- SurAll = Region[{GammaD, GammaN, GammaR, GammaBT}] ;
- VolAll = Region[{Omega, OmegaPml}] ;
- TotAll = Region[{VolAll, SurAll}] ;
+ DefineGroup[ Omega, OmegaPml, GammaD, GammaN, GammaR, GammaBT ];
+ SurAll = Region[{GammaD, GammaN, GammaR, GammaBT}];
+ VolAll = Region[{Omega, OmegaPml}];
+ TotAll = Region[{VolAll, SurAll}];
}
-Jacobian {
- { Name JVol ;
- Case {
- { Region All ; Jacobian Vol ; }
- }
- }
- { Name JSur ;
- Case {
- { Region All ; Jacobian Sur ; }
- }
- }
-}
-
-Integration {
- { Name I1 ;
- Case {
- { Type Gauss ;
- Case {
- { GeoElement Point ; NumberOfPoints 1 ; }
- { GeoElement Line ; NumberOfPoints 4 ; }
- { GeoElement Triangle ; NumberOfPoints 6 ; }
- { GeoElement Quadrangle ; NumberOfPoints 7 ; }
- { GeoElement Tetrahedron ; NumberOfPoints 15 ; }
- { GeoElement Hexahedron ; NumberOfPoints 34 ; }
- }
- }
- }
- }
-}
+//==============================================================================
FunctionSpace {
- { Name pSpace ; Type Form0 ;
+ { Name pSpace; Type Form0;
BasisFunction {
- { Name sn ; NameOfCoef un ; Function BF_Node ;
- Support TotAll ; Entity NodesOf[All] ; }
- If (ORDER==2)
- { Name sn2 ; NameOfCoef un2 ; Function BF_Node_2E ;
- Support TotAll ; Entity EdgesOf[All] ; }
- EndIf
+ { Name sn; NameOfCoef un; Function BF_Node;
+ Support TotAll; Entity NodesOf[All]; }
}
Constraint {
- { NameOfCoef un ; EntityType NodesOf ; NameOfConstraint pConstraint ; }
- If (ORDER==2)
- { NameOfCoef un2 ; EntityType EdgesOf ; NameOfConstraint pConstraint2 ; }
- EndIf
+ { NameOfCoef un; EntityType NodesOf; NameOfConstraint pConstraint; }
}
}
- { Name vSpace ; Type Form0 ;
+ { Name vSpace; Type Form0;
BasisFunction {
- { Name sn ; NameOfCoef vnA ; Function BF_Node ;
- Support TotAll ; Entity NodesOf[All] ; }
+ { Name sn; NameOfCoef vnA; Function BF_Node;
+ Support TotAll; Entity NodesOf[All]; }
}
Constraint {
- { NameOfCoef vnA ; EntityType NodesOf ; NameOfConstraint uConstraint ; }
+ { NameOfCoef vnA; EntityType NodesOf; NameOfConstraint uConstraint; }
}
}
- If (DIM==2)
- { Name uSpace ; Type Form2P ;
+If (FLAG_DIM==2)
+ { Name uSpace; Type Form2P;
BasisFunction {
- { Name se ; NameOfCoef he ; Function BF_PerpendicularFacet ;
- Support TotAll ; Entity EdgesOf[All] ; }
+ { Name se; NameOfCoef he; Function BF_PerpendicularFacet;
+ Support TotAll; Entity EdgesOf[All]; }
}
Constraint {
- { NameOfCoef he ; EntityType EdgesOf ; NameOfConstraint uConstraint ; }
+ { NameOfCoef he; EntityType EdgesOf; NameOfConstraint uConstraint; }
}
}
- EndIf
- If (DIM==3)
- { Name uSpace ; Type Form2 ;
+EndIf
+If (FLAG_DIM==3)
+ { Name uSpace; Type Form2;
BasisFunction {
- { Name se ; NameOfCoef he ; Function BF_Facet ;
- Support TotAll ; Entity FacetsOf[All] ; }
+ { Name se; NameOfCoef he; Function BF_Facet;
+ Support TotAll; Entity FacetsOf[All]; }
}
Constraint {
- { NameOfCoef he ; EntityType FacetsOf ; NameOfConstraint uConstraint ; }
+ { NameOfCoef he; EntityType FacetsOf; NameOfConstraint uConstraint; }
}
}
- EndIf
+EndIf
}
+//==============================================================================
+
Formulation {
- { Name WaveEquation_Form ; Type FemEquation ;
+ { Name WaveEqn_Form; Type FemEquation;
Quantity {
- { Name p ; Type Local ; NameOfSpace pSpace ; }
+ { Name p; Type Local; NameOfSpace pSpace; }
}
Equation {
- Galerkin { DtDtDof [ Dof{p} , {p} ] ;
- In VolAll ; Integration I1 ; Jacobian JVol ; }
- Galerkin { [ Dof{d p} , {d p} ] ;
- In Omega ; Jacobian JVol ; Integration I1 ; }
+ Galerkin { DtDtDof [ Dof{p} , {p} ];
+ In Omega; Jacobian JVol; Integration I1; }
+ Galerkin { [ Dof{d p} , {d p} ];
+ In Omega; Jacobian JVol; Integration I1; }
// BC: Neumann
- Galerkin { [ -fNeumann[] , {p} ] ;
- In GammaN ; Jacobian JSur ; Integration I1 ; }
+ Galerkin { [ -fNeumann[] , {p} ];
+ In GammaN; Jacobian JSur; Integration I1; }
// BC: Radiation
- // Galerkin { [ -I[]*k * Dof{p} , {p} ] ;
- Galerkin { DtDof [ -Dof{p} , {p} ] ;
- In GammaR ; Jacobian JSur ; Integration I1 ; }
+ Galerkin { DtDof [ Dof{p} , {p} ];
+ In GammaR; Jacobian JSur; Integration I1; }
// BC: Bayliss-Turkel
- // Galerkin { [ -I[]*k * Dof{p} , {p} ] ;
- Galerkin { DtDof [ Dof{p} , {p} ] ;
- In GammaBT ; Jacobian JSur ; Integration I1 ; }
- Galerkin { [ alphaBT[] * Dof{p} , {p} ] ;
- In GammaBT ; Jacobian JSur ; Integration I1 ; }
- Galerkin { [ betaBT[] * Dof{d p} , {d p} ] ;
- In GammaBT ; Jacobian JSur ; Integration I1 ; }
+ Galerkin { DtDof [ Dof{p} , {p} ];
+ In GammaBT; Jacobian JSur; Integration I1; }
+ Galerkin { [ alphaBT[] * Dof{p} , {p} ];
+ In GammaBT; Jacobian JSur; Integration I1; }
+ Galerkin { [ betaBT[] * Dof{d p} , {d p} ];
+ In GammaBT; Jacobian JSur; Integration I1; }
// PML (only for resolutions in the frequency domain)
- Galerkin { [ pmlTens[] * Dof{d p} , {d p} ] ;
- In OmegaPml ; Jacobian JVol ; Integration I1 ; }
- Galerkin { [ -k^2 * (1/pmlScal[]) * Dof{p} , {p} ] ;
- In OmegaPml ; Jacobian JVol ; Integration I1 ; }
+ Galerkin { [ pmlTens[] * Dof{d p} , {d p} ];
+ In OmegaPml; Jacobian JVol; Integration I1; }
+ Galerkin { [ -k^2 * (1/pmlScal[]) * Dof{p} , {p} ];
+ In OmegaPml; Jacobian JVol; Integration I1; }
}
}
- { Name MembraneEquations_Form ; Type FemEquation ;
+ { Name MembraneSys_Form; Type FemEquation;
Quantity {
- { Name p ; Type Local ; NameOfSpace pSpace ; }
- { Name v ; Type Local ; NameOfSpace vSpace ; }
+ { Name p; Type Local; NameOfSpace pSpace; }
+ { Name v; Type Local; NameOfSpace vSpace; }
}
Equation {
- Galerkin { DtDof [ Dof{p} , {p} ] ;
- In VolAll ; Integration I1 ; Jacobian JVol ; }
- Galerkin { [ - Dof{v} , {p} ] ;
- In VolAll ; Integration I1 ; Jacobian JVol ; }
+ Galerkin { DtDof [ Dof{p} , {p} ];
+ In VolAll; Integration I1; Jacobian JVol; }
+ Galerkin { [ -Dof{v} , {p} ];
+ In VolAll; Integration I1; Jacobian JVol; }
- Galerkin { DtDof [ Dof{v} , {v} ] ;
- In VolAll ; Integration I1 ; Jacobian JVol ; }
- Galerkin { [ Dof{d p} , {d v} ] ;
- In VolAll ; Integration I1 ; Jacobian JVol ; }
+ Galerkin { DtDof [ Dof{v} , {v} ];
+ In VolAll; Integration I1; Jacobian JVol; }
+ Galerkin { [ Dof{d p} , {d v} ];
+ In VolAll; Integration I1; Jacobian JVol; }
}
}
- { Name EulerEquations_Form ; Type FemEquation ;
+ { Name EulerSys_Form; Type FemEquation;
Quantity {
- { Name p ; Type Local ; NameOfSpace pSpace ; }
- { Name u ; Type Local ; NameOfSpace uSpace ; }
+ { Name p; Type Local; NameOfSpace pSpace; }
+ { Name u; Type Local; NameOfSpace uSpace; }
}
Equation {
- Galerkin { DtDof [ Dof{p} , {p} ] ;
- In VolAll ; Integration I1 ; Jacobian JVol ; }
- Galerkin { [ - Dof{u} , {d p} ] ;
- In VolAll ; Integration I1 ; Jacobian JVol ; }
+ Galerkin { DtDof [ Dof{p} , {p} ];
+ In VolAll; Integration I1; Jacobian JVol; }
+ Galerkin { [ - Dof{u} , {d p} ];
+ In VolAll; Integration I1; Jacobian JVol; }
- Galerkin { DtDof [ Dof{u} , {u} ] ;
- In VolAll ; Integration I1 ; Jacobian JVol ; }
- Galerkin { [ - Dof{p} , {d u} ] ;
- In VolAll ; Integration I1 ; Jacobian JVol ; }
+ Galerkin { DtDof [ Dof{u} , {u} ];
+ In VolAll; Integration I1; Jacobian JVol; }
+ Galerkin { [ - Dof{p} , {d u} ];
+ In VolAll; Integration I1; Jacobian JVol; }
// BC: Radiation
- Galerkin { [ Dof{p} , {p} ] ;
- In GammaR ; Integration I1 ; Jacobian JSur ; }
- Galerkin { [ (Dof{u} * Normal[]) * Normal[] , {u} ] ;
- In GammaR ; Integration I1 ; Jacobian JSur ; }
-
+ Galerkin { [ Dof{p} , {p} ];
+ In GammaR; Integration I1; Jacobian JSur; }
+ Galerkin { [ (Dof{u} * Normal[]) * Normal[] , {u} ];
+ In GammaR; Integration I1; Jacobian JSur; }
}
}
}
+//==============================================================================
+
Resolution {
- { Name WaveEquation_FreqReso ;
+If(FLAG_RES == RES_FREQ)
+ { Name WaveEqn_Reso;
System {
- { Name A ; NameOfFormulation WaveEquation_Form ; Frequency Freq ; }
+ { Name A; NameOfFormulation WaveEqn_Form; Frequency FREQ; }
}
Operation {
- CreateDir["res/"] ;
- Generate[A] ;
- Solve[A] ;
- SaveSolution[A] ;
+ CreateDir["res/"];
+ Generate[A]; Solve[A]; SaveSolution[A];
}
}
- { Name MembraneEquations_FreqReso ;
+ { Name MembraneSys_Reso;
System {
- { Name A ; NameOfFormulation MembraneEquations_Form ; Frequency Freq ; }
+ { Name A; NameOfFormulation MembraneSys_Form; Frequency FREQ; }
}
Operation {
- CreateDir["res/"] ;
- Generate[A] ;
- Solve[A] ;
- SaveSolution[A] ;
+ CreateDir["res/"];
+ Generate[A]; Solve[A]; SaveSolution[A];
}
}
- { Name EulerEquations_FreqReso ;
+ { Name EulerSys_Reso;
System {
- { Name A ; NameOfFormulation EulerEquations_Form ; Frequency Freq ; }
+ { Name A; NameOfFormulation EulerSys_Form; Frequency FREQ; }
}
Operation {
- CreateDir["res/"] ;
- Generate[A] ;
- Solve[A] ;
- SaveSolution[A] ;
+ CreateDir["res/"];
+ Generate[A]; Solve[A]; SaveSolution[A];
}
}
- { Name WaveEquation_TimeReso ;
+EndIf
+If(FLAG_RES == RES_TIME)
+ { Name WaveEqn_Reso;
System {
- { Name A ; NameOfFormulation WaveEquation_Form ; }
+ { Name A; NameOfFormulation WaveEqn_Form; }
}
Operation {
- CreateDir["res/"] ;
- InitSolution[A] ;
- InitSolution[A] ;
+ CreateDir["res/"];
+ InitSolution[A];
+ InitSolution[A];
TimeLoopNewmark[0., tf, dt, beta, gamma] {
- Generate[A] ;
- Solve[A] ;
- SaveSolution[A] ;
- }
+ Generate[A]; Solve[A]; SaveSolution[A]; }
}
}
- { Name MembraneEquations_TimeReso ;
+ { Name MembraneSys_Reso;
System {
- { Name A ; NameOfFormulation MembraneEquations_Form ; }
+ { Name A; NameOfFormulation MembraneSys_Form; }
}
Operation {
- CreateDir["res/"] ;
- InitSolution[A] ;
+ CreateDir["res/"];
+ InitSolution[A];
TimeLoopTheta{
- Time0 0 ; DTime dt ; Theta theta ; TimeMax tf ;
- Operation {
- Generate[A] ;
- Solve[A] ;
- SaveSolution[A] ;
- }
+ Time0 0; DTime dt; Theta theta; TimeMax tf;
+ Operation { Generate[A]; Solve[A]; SaveSolution[A]; }
}
}
}
- { Name EulerEquations_TimeReso ;
+ { Name EulerSys_Reso;
System {
- { Name A ; NameOfFormulation EulerEquations_Form ; }
+ { Name A; NameOfFormulation EulerSys_Form; }
}
Operation {
- CreateDir["res/"] ;
- InitSolution[A] ;
+ CreateDir["res/"];
+ InitSolution[A];
TimeLoopTheta{
- Time0 0 ; DTime dt ; Theta theta ; TimeMax tf ;
- Operation {
- Generate[A] ;
- Solve[A] ;
- SaveSolution[A] ;
- }
+ Time0 0; DTime dt; Theta theta; TimeMax tf;
+ Operation { Generate[A]; Solve[A]; SaveSolution[A]; }
}
}
}
+EndIf
+}
+
+//==============================================================================
+
+PostProcessing {
+ { Name WaveEqn_PostPro; NameOfFormulation WaveEqn_Form; NameOfSystem A;
+ Quantity {
+ { Name p; Value{ Local{ [ {p} ]; In VolAll; Jacobian JVol; }}}
+ }
+ }
+ { Name MembraneSys_PostPro; NameOfFormulation MembraneSys_Form; NameOfSystem A;
+ Quantity {
+ { Name p; Value{ Local{ [ {p} ]; In VolAll; Jacobian JVol; }}}
+ { Name v; Value{ Local{ [ {v} ]; In VolAll; Jacobian JVol; }}}
+ }
+ }
+ { Name EulerSys_PostPro; NameOfFormulation EulerSys_Form; NameOfSystem A;
+ Quantity {
+ { Name p; Value{ Local{ [ {p} ]; In VolAll; Jacobian JVol; }}}
+ { Name u; Value{ Local{ [ {u} ]; In VolAll; Jacobian JVol; }}}
+ }
+ }
+}
+
+//==============================================================================
+
+PostOperation {
+ { Name WaveEqn_PostOp; NameOfPostProcessing WaveEqn_PostPro;
+ Operation {
+ Print[ p, OnElementsOf VolAll, File "res/p.pos"];
+ }
+ }
+ { Name MembraneSys_PostOp; NameOfPostProcessing MembraneSys_PostPro;
+ Operation {
+ Print[ p, OnElementsOf VolAll, File "res/p.pos"];
+ Print[ v, OnElementsOf VolAll, File "res/u.pos"];
+ }
+ }
+ { Name EulerSys_PostOp; NameOfPostProcessing EulerSys_PostPro;
+ Operation {
+ Print[ p, OnElementsOf VolAll, File "res/p.pos"];
+ Print[ u, OnElementsOf VolAll, File "res/u.pos"];
+ }
+ }
}
diff --git a/AcademicWaves/formulations_vectorWaves.pro b/AcademicWaves/formulations_vectorWaves.pro
index ce4c76d..1c703b7 100644
--- a/AcademicWaves/formulations_vectorWaves.pro
+++ b/AcademicWaves/formulations_vectorWaves.pro
@@ -5,51 +5,21 @@
//========================================================
Function{
- DefineFunction[ fNeumann, fDirichlet ] ;
+ DefineFunction[ fNeumann, fDirichlet ];
}
Group{
- DefineGroup[ Omega, GammaD, GammaN, GammaR ] ;
- SurAll = Region[{GammaD, GammaN, GammaR }] ;
- VolAll = Region[{Omega}] ;
- TotAll = Region[{VolAll,SurAll}] ;
-}
-
-Jacobian {
- { Name JVol ;
- Case {
- { Region All ; Jacobian Vol ; }
- }
- }
- { Name JSur ;
- Case {
- { Region All ; Jacobian Sur ; }
- }
- }
-}
-
-Integration {
- { Name I1 ;
- Case {
- { Type Gauss ;
- Case {
- { GeoElement Line ; NumberOfPoints 4 ; }
- { GeoElement Triangle ; NumberOfPoints 4 ; }
- { GeoElement Quadrangle ; NumberOfPoints 4 ; }
- { GeoElement Tetrahedron ; NumberOfPoints 4 ; }
- { GeoElement Hexahedron ; NumberOfPoints 6 ; }
- { GeoElement Prism ; NumberOfPoints 9 ; }
- }
- }
- }
- }
+ DefineGroup[ Omega, GammaD, GammaN, GammaR ];
+ SurAll = Region[{GammaD, GammaN, GammaR }];
+ VolAll = Region[{Omega}];
+ TotAll = Region[{VolAll,SurAll}];
}
FunctionSpace {
- { Name eSpace ; Type Form1 ;
+ { Name eSpace; Type Form1;
BasisFunction {
- { Name se ; NameOfCoef ee ; Function BF_Edge ;
- Support TotAll ; Entity EdgesOf[All] ; }
+ { Name se; NameOfCoef ee; Function BF_Edge;
+ Support TotAll; Entity EdgesOf[All]; }
}
}
{ Name eSpace_DirichletBC; Type Form1;
@@ -61,16 +31,16 @@ FunctionSpace {
}
Formulation {
- { Name Form ; Type FemEquation ;
+ { Name Form; Type FemEquation;
Quantity {
- { Name e ; Type Local ; NameOfSpace eSpace ; }
- { Name lambda ; Type Local ; NameOfSpace eSpace_DirichletBC ; }
+ { Name e; Type Local; NameOfSpace eSpace; }
+ { Name lambda; Type Local; NameOfSpace eSpace_DirichletBC; }
}
Equation {
- Galerkin { DtDtDof [ epsilon[] * Dof{e} , {e} ] ;
- In VolAll ; Integration I1 ; Jacobian JVol ; }
- Galerkin { [ nu[] * Dof{d e} , {d e} ] ;
- In VolAll ; Integration I1 ; Jacobian JVol ; }
+ Galerkin { DtDtDof [ epsilon[] * Dof{e} , {e} ];
+ In VolAll; Integration I1; Jacobian JVol; }
+ Galerkin { [ nu[] * Dof{d e} , {d e} ];
+ In VolAll; Integration I1; Jacobian JVol; }
// BC: Dirichlet
Galerkin { [ Dof{lambda} , {e} ];
@@ -81,21 +51,21 @@ Formulation {
In GammaD; Integration I1; Jacobian JSur; }
// BC: Radiation
- Galerkin { [ I[]*k * nu[] * ( extNormal[] /\ Dof{e} ) /\ extNormal[] , {e} ] ;
- In GammaR ; Integration I1 ; Jacobian JSur ; }
+ Galerkin { [ I[]*k * nu[] * ( extNormal[] /\ Dof{e} ) /\ extNormal[] , {e} ];
+ In GammaR; Integration I1; Jacobian JSur; }
}
}
}
Resolution {
- { Name Reso ;
+ { Name Reso;
System {
- { Name A ; NameOfFormulation Form ; Type Complex ; Frequency Freq ; }
+ { Name A; NameOfFormulation Form; Type Complex; Frequency FREQ; }
}
Operation {
- Generate[A] ;
- Solve[A] ;
- SaveSolution[A] ;
+ Generate[A];
+ Solve[A];
+ SaveSolution[A];
}
}
}
diff --git a/AcademicWaves/gaussInCuboid.dat b/AcademicWaves/gaussInCuboid.dat
index 6ff388d..1f665cc 100644
--- a/AcademicWaves/gaussInCuboid.dat
+++ b/AcademicWaves/gaussInCuboid.dat
@@ -3,6 +3,9 @@
// File: parameters
//========================================================
+FLAG_DIM = 3;
+FLAG_RES = RES_TIME;
+
//---------------------
// Physical parameters
//---------------------
@@ -15,16 +18,12 @@ tf = 1;
// Numerical parameters
//---------------------
-lc = 0.1;
-Flag_SecondOrder = 0;
-dt = lc/4;
-
-// For the theta-scheme (if wave system)
-theta = 0.5;
+LC = 0.1;
+dt = LC/4;
-// For the Newmark scheme (if wave equation)
-beta = 0.25;
-gamma = 0.5;
+theta = 0.5; // for theta-scheme (if wave system)
+beta = 0.25; // for Newmark scheme (if wave equation)
+gamma = 0.5; // for Newmark scheme (if wave equation)
//---------------------
// Geometrical tags
diff --git a/AcademicWaves/gaussInCuboid.geo b/AcademicWaves/gaussInCuboid.geo
index 97cc4aa..d424cdf 100644
--- a/AcademicWaves/gaussInCuboid.geo
+++ b/AcademicWaves/gaussInCuboid.geo
@@ -5,7 +5,9 @@
Include "gaussInCuboid.dat";
-Mesh.CharacteristicLengthMax = lc;
+Mesh.CharacteristicLengthMax = LC;
+Mesh.CharacteristicLengthFactor = 1;
+Mesh.Optimize = 1;
// Points
diff --git a/AcademicWaves/gaussInCuboid.pro b/AcademicWaves/gaussInCuboid.pro
index 89a0de1..0609882 100644
--- a/AcademicWaves/gaussInCuboid.pro
+++ b/AcademicWaves/gaussInCuboid.pro
@@ -3,49 +3,39 @@
// File: GetDP simulation
//========================================================
-Include "gaussInCuboid.dat" ;
+Include "gaussInCuboid.dat";
Group {
- GammaD = Region[{}] ;
- GammaS = Region[{BOUNDARY}] ;
- Omega = Region[{DOMAIN}] ;
+ GammaD = Region[{}];
+ GammaN = Region[{}];
+ GammaR = Region[{BOUNDARY}];
+ Omega = Region[{DOMAIN}];
}
Function {
- InitValue[] = Exp[- ((X[])^2 + (Y[])^2 + (Z[])^2) / (r*r)] ;
+ InitValue[] = Exp[-(X[]^2+Y[]^2+Z[]^2)/(r*r)];
}
Constraint {
- { Name pConstraint ;
+ { Name pConstraint;
Case {
- { Region Omega ; Type Init ; Value InitValue[] ; }
+ { Region Omega; Type Init; Value InitValue[]; }
}
}
- { Name uConstraint ;
+ { Name uConstraint;
Case {
- { Region GammaD ; Type Assign ; Value 0. ; }
- { Region Omega ; Type Init ; Value 0. ; }
+ { Region GammaD; Type Assign; Value 0.; }
+ { Region Omega; Type Init; Value 0.; }
}
}
}
-Include "form_scalarWaveSyst_time3D.pro" ;
-//Include "form_scalarWaveEqn_time3D.pro" ;
+Include "formulations_scalarWaves.pro";
-PostProcessing {
- { Name PostPro ; NameOfFormulation Form ; NameOfSystem A ;
- Quantity {
- { Name p ; Value{ Local{ [ {p} ] ; In VolAll ; } } }
- { Name u ; Value{ Local{ [ {u} ] ; In VolAll ; } } }
- }
- }
-}
-
-PostOperation {
- { Name PostOp ; NameOfPostProcessing PostPro ;
- Operation {
- Print[ p, OnElementsOf VolAll, File "res/p.pos"] ;
- Print[ u, OnElementsOf VolAll, File "res/u.pos"] ;
- }
- }
-}
+DefineConstant[
+ C_ = {"-solve -pos -v2", Name "GetDP/9ComputeCommand", Visible 0},
+ R_ = {"WaveEqn_Reso", Name "GetDP/1ResolutionChoices", Visible 1,
+ Choices {"WaveEqn_Reso", "MembraneSys_Reso", "EulerSys_Reso"}},
+ P_ = {"WaveEqn_PostOp", Name "GetDP/2PostOperationChoices", Visible 1,
+ Choices {"WaveEqn_PostOp", "MembraneSys_PostOp", "EulerSys_PostOp"}}
+];
diff --git a/AcademicWaves/gaussInSquare.dat b/AcademicWaves/gaussInSquare.dat
index daa4ecc..c669e0b 100644
--- a/AcademicWaves/gaussInSquare.dat
+++ b/AcademicWaves/gaussInSquare.dat
@@ -3,28 +3,27 @@
// File: parameters
//========================================================
+FLAG_DIM = 2;
+FLAG_RES = RES_TIME;
+
//---------------------
// Physical parameters
//---------------------
ldom = 1;
r = ldom/5 ;
-tf = 1;
+tf = 1.;
//---------------------
// Numerical parameters
//---------------------
-lc = 0.025;
-Flag_SecondOrder = 0;
-dt = lc/4;
-
-// For the theta-scheme (if wave system)
-theta = 0.5;
+LC = 0.025;
+dt = LC/4;
-// For the Newmark scheme (if wave equation)
-beta = 0.25;
-gamma = 0.5;
+theta = 0.5; // for theta-scheme (if wave system)
+beta = 0.25; // for Newmark scheme (if wave equation)
+gamma = 0.5; // for Newmark scheme (if wave equation)
//---------------------
// Geometrical tags
diff --git a/AcademicWaves/gaussInSquare.geo b/AcademicWaves/gaussInSquare.geo
index fd89fe5..9d9b35c 100644
--- a/AcademicWaves/gaussInSquare.geo
+++ b/AcademicWaves/gaussInSquare.geo
@@ -5,10 +5,14 @@
Include "gaussInSquare.dat";
-Point(1) = {-ldom/2, -ldom/2, 0, lc};
-Point(2) = { ldom/2, -ldom/2, 0, lc};
-Point(3) = { ldom/2, ldom/2, 0, lc};
-Point(4) = {-ldom/2, ldom/2, 0, lc};
+Mesh.CharacteristicLengthMax = LC;
+Mesh.CharacteristicLengthFactor = 1;
+Mesh.Optimize = 1;
+
+Point(1) = {-ldom/2, -ldom/2, 0};
+Point(2) = { ldom/2, -ldom/2, 0};
+Point(3) = { ldom/2, ldom/2, 0};
+Point(4) = {-ldom/2, ldom/2, 0};
Line(1) = {1, 2};
Line(2) = {2, 3};
diff --git a/AcademicWaves/gaussInSquare.pro b/AcademicWaves/gaussInSquare.pro
index 8c11a0a..6629829 100644
--- a/AcademicWaves/gaussInSquare.pro
+++ b/AcademicWaves/gaussInSquare.pro
@@ -3,52 +3,39 @@
// File: GetDP simulation
//========================================================
-Include "gaussInSquare.dat" ;
+Include "gaussInSquare.dat";
Group {
- GammaD = Region[{}] ;
- GammaN = Region[{BOUNDARY}] ;
- GammaS = Region[{}] ;
- Omega = Region[{DOMAIN}] ;
+ GammaD = Region[{}];
+ GammaN = Region[{}];
+ GammaR = Region[{BOUNDARY}];
+ Omega = Region[{DOMAIN}];
}
Function {
- InitValue[] = Exp[- ((X[])^2 + (Y[])^2) / (r*r)] ;
- fNeumann[] = 0.;
+ InitValue[] = Exp[-(X[]^2+Y[]^2)/(r*r)];
}
Constraint {
- { Name pConstraint ;
+ { Name pConstraint;
Case {
- { Region Omega ; Type Init ; Value InitValue[] ; }
+ { Region Omega; Type Init; Value InitValue[]; }
}
}
- { Name uConstraint ;
+ { Name uConstraint;
Case {
- { Region GammaD ; Type Assign ; Value 0. ; }
- { Region Omega ; Type Init ; Value 0. ; }
+ { Region GammaD; Type Assign; Value 0.; }
+ { Region Omega; Type Init; Value 0.; }
}
}
}
-//Include "form_scalarWaveSyst_time2D.pro" ;
-//Include "form_scalarWaveSyst_time2D_membrane.pro" ;
-Include "form_scalarWaveEqn_time.pro" ;
+Include "formulations_scalarWaves.pro";
-PostProcessing {
- { Name PostPro ; NameOfFormulation Form ; NameOfSystem A ;
- Quantity {
- { Name p ; Value{ Local{ [ {p} ] ; In VolAll ; } } }
- { Name u ; Value{ Local{ [ {u} ] ; In VolAll ; } } }
- }
- }
-}
-
-PostOperation {
- { Name PostOp ; NameOfPostProcessing PostPro ;
- Operation {
- Print[ p, OnElementsOf VolAll, File "res/p.pos"] ;
- Print[ u, OnElementsOf VolAll, File "res/u.pos"] ;
- }
- }
-}
+DefineConstant[
+ C_ = {"-solve -pos -v2", Name "GetDP/9ComputeCommand", Visible 0},
+ R_ = {"WaveEqn_Reso", Name "GetDP/1ResolutionChoices", Visible 1,
+ Choices {"WaveEqn_Reso", "MembraneSys_Reso", "EulerSys_Reso"}},
+ P_ = {"WaveEqn_PostOp", Name "GetDP/2PostOperationChoices", Visible 1,
+ Choices {"WaveEqn_PostOp", "MembraneSys_PostOp", "EulerSys_PostOp"}}
+];
diff --git a/AcademicWaves/loudspeaker.dat b/AcademicWaves/loudspeaker.dat
index 302f554..04de0cf 100644
--- a/AcademicWaves/loudspeaker.dat
+++ b/AcademicWaves/loudspeaker.dat
@@ -3,23 +3,13 @@
// File: parameters
//========================================================
+FLAG_DIM = 2;
+FLAG_RES = RES_FREQ;
//---------------------
-// Model choices
+// Geometrical parameters
//---------------------
-DefineConstant[
- Flag_ModelEquations = {1, Choices{1="Wave equation (p)", 2="Wave system (p and u)"},
- Name "Input/01Model equations", Highlight "Black"}
-];
-
-
-//---------------------
-// Physical parameters
-//---------------------
-
-// Geometry
-
feed = 0.1;
app = 0.5;
a = 0.1;
@@ -28,24 +18,15 @@ thick = 0.03;
x0 = thick-(a+b+thick)/2;
Ro = 1;
+//---------------------
// Incident signal
+//---------------------
k = 20;
-lambda = 2*Pi/k;
-nlambda = 10;
+FREQ = k/(2*Pi);
+LC = 0.01;
theta = 0;
-
-//---------------------
-// Numerical parameters
-//---------------------
-
-lc0 = 0.02;
-lc1 = 0.02;
-
-ORDER = 1;
-
-
//---------------------
// Geometrical tags
//---------------------
diff --git a/AcademicWaves/loudspeaker.geo b/AcademicWaves/loudspeaker.geo
index a8e3a17..6dbd37d 100644
--- a/AcademicWaves/loudspeaker.geo
+++ b/AcademicWaves/loudspeaker.geo
@@ -5,24 +5,27 @@
Include "loudspeaker.dat";
+Mesh.CharacteristicLengthMax = LC;
+Mesh.CharacteristicLengthFactor = 1;
+Mesh.Optimize = 1;
//--------------------------------------------------------
// The loudspeaker
//--------------------------------------------------------
-cen0=newp ; Point(newp) = {0,0,0,lc0} ;
-pObs[]+=newp ; Point(newp) = {x0, -feed/2, 0, lc0} ;
-pObs[]+=newp ; Point(newp) = {x0, feed/2, 0, lc0} ;
-pObs[]+=newp ; Point(newp) = {a+x0, feed/2, 0, lc0} ;
-pObs[]+=newp ; Point(newp) = {a+b+x0, app/2, 0, lc0} ;
-pObs[]+=newp ; Point(newp) = {a+b+x0, app/2+thick, 0, lc0} ;
-pObs[]+=newp ; Point(newp) = {a+x0, feed/2+thick, 0, lc0} ;
-pObs[]+=newp ; Point(newp) = {-thick+x0, feed/2+thick, 0, lc0} ;
-pObs[]+=newp ; Point(newp) = {-thick+x0, -feed/2-thick, 0, lc0} ;
-pObs[]+=newp ; Point(newp) = {a+x0, -feed/2-thick, 0, lc0} ;
-pObs[]+=newp ; Point(newp) = {a+b+x0, -app/2-thick, 0, lc0} ;
-pObs[]+=newp ; Point(newp) = {a+b+x0, -app/2, 0, lc0} ;
-pObs[]+=newp ; Point(newp) = {a+x0, -feed/2, 0, lc0} ;
+cen0=newp ; Point(newp) = {0,0,0} ;
+pObs[]+=newp ; Point(newp) = {x0, -feed/2, 0} ;
+pObs[]+=newp ; Point(newp) = {x0, feed/2, 0} ;
+pObs[]+=newp ; Point(newp) = {a+x0, feed/2, 0} ;
+pObs[]+=newp ; Point(newp) = {a+b+x0, app/2, 0} ;
+pObs[]+=newp ; Point(newp) = {a+b+x0, app/2+thick, 0} ;
+pObs[]+=newp ; Point(newp) = {a+x0, feed/2+thick, 0} ;
+pObs[]+=newp ; Point(newp) = {-thick+x0, feed/2+thick, 0} ;
+pObs[]+=newp ; Point(newp) = {-thick+x0, -feed/2-thick, 0} ;
+pObs[]+=newp ; Point(newp) = {a+x0, -feed/2-thick, 0} ;
+pObs[]+=newp ; Point(newp) = {a+b+x0, -app/2-thick, 0} ;
+pObs[]+=newp ; Point(newp) = {a+b+x0, -app/2, 0} ;
+pObs[]+=newp ; Point(newp) = {a+x0, -feed/2, 0} ;
For i In {0:#pObs[]-1}
lObs[]+=newl ; Line(newl) = {pObs[i],pObs[(i+1)%#pObs[]]} ;
EndFor
@@ -32,16 +35,15 @@ llObs=newll ; Line Loop(newll) = lObs[] ;
Physical Line(GAMMA_SOURCE) = {-lObs[{0}]} ;
Physical Line(GAMMA_WALL) = {-lObs[{1:11}]} ;
-
//--------------------------------------------------------
// Main domain with its boundary
//--------------------------------------------------------
-cen0 = newp ; Point(newp) = {0,0,0,lc0} ;
-pi[]+=newp ; Point(newp) = {Ro,0,0,lc1} ;
-pi[]+=newp ; Point(newp) = {0,Ro,0,lc1} ;
-pi[]+=newp ; Point(newp) = {-Ro,0,0,lc1} ;
-pi[]+=newp ; Point(newp) = {0,-Ro,0,lc1} ;
+cen0 = newp ; Point(newp) = {0,0,0} ;
+pi[]+=newp ; Point(newp) = {Ro,0,0} ;
+pi[]+=newp ; Point(newp) = {0,Ro,0} ;
+pi[]+=newp ; Point(newp) = {-Ro,0,0} ;
+pi[]+=newp ; Point(newp) = {0,-Ro,0} ;
li[]+=newl ; Circle(newl) = {pi[0], cen0, pi[1]} ;
li[]+=newl ; Circle(newl) = {pi[1], cen0, pi[2]} ;
li[]+=newl ; Circle(newl) = {pi[2], cen0, pi[3]} ;
diff --git a/AcademicWaves/loudspeaker.pro b/AcademicWaves/loudspeaker.pro
index 0322f9a..85346e9 100644
--- a/AcademicWaves/loudspeaker.pro
+++ b/AcademicWaves/loudspeaker.pro
@@ -3,92 +3,39 @@
// File: GetDP simulation
//========================================================
-Include "loudspeaker.dat" ;
+Include "loudspeaker.dat";
Group {
- Omega = Region[{DOMAIN}] ;
- If (Flag_ModelEquations==1)
- GammaD = Region[{GAMMA_SOURCE}] ;
- GammaN = Region[{GAMMA_WALL}] ;
- EndIf
- If (Flag_ModelEquations==2)
- GammaD_p = Region[{GAMMA_SOURCE}] ;
- GammaD_u = Region[{GAMMA_WALL}] ;
- GammaD = Region[{GammaD_p,GammaD_u}] ;
- EndIf
- GammaR = Region[{GAMMA_DOM}] ;
+ Omega = Region[{DOMAIN}];
+ GammaD = Region[{GAMMA_SOURCE}];
+ GammaN = Region[{GAMMA_WALL}];
+ GammaR = Region[{GAMMA_DOM}];
}
Function {
- I[] = Complex[0, 1] ;
+ I[] = Complex[0, 1];
// Incident signal
- kVec[] = Vector[ k*Cos[theta], k*Sin[theta], 0] ;
- pInc[] = Complex[ Cos[kVec[]*XYZ[]], Sin[kVec[]*XYZ[]] ] ;
+ kVec[] = Vector[ k*Cos[theta], k*Sin[theta], 0];
+ pInc[] = Complex[ Cos[kVec[]*XYZ[]], Sin[kVec[]*XYZ[]] ];
// BC: Neumann
- If (Flag_ModelEquations==2)
- fNeumann[] = 0. ;
- EndIf
+ fNeumann[] = 0.;
}
Constraint {
- { Name pConstraint ;
- Case {
- { Region Region[{GAMMA_SOURCE}] ; Value pInc[] ; }
- }
+ { Name pConstraint;
+ Case {{ Region GammaD; Value pInc[]; }}
}
- If(Flag_SecondOrder)
- { Name pConstraint2 ;
- Case {
- { Region Region[{GAMMA_SOURCE}] ; Value 0. ; }
- }
+ { Name uConstraint;
+ Case {{ Region GammaN; Value 0.; }}
}
- EndIf
- { Name uConstraint ;
- Case {
- { Region Region[{GAMMA_WALL}] ; Value 0. ; }
- }
- }
- If(Flag_SecondOrder)
- { Name uConstraint2;
- Case {
- { Region Region[{GAMMA_WALL}] ; Value 0. ; }
- }
- }
- EndIf
}
-Include "formulations_scalarWaves.pro" ;
-
-PostProcessing {
- If (Flag_ModelEquations==1)
- { Name PostPro ; NameOfFormulation WaveEquation_Form ;
- Quantity {
- { Name p ; Value{ Local{ [ {p} ] ; In VolAll ; Jacobian JVol ; } } }
- }
- }
- EndIf
- If (Flag_ModelEquations==2)
- { Name PostPro ; NameOfFormulation EulerEquations_Form ;
- Quantity {
- { Name p ; Value{ Local{ [ {p} ] ; In VolAll ; Jacobian JVol ; } } }
- }
- }
- EndIf
-}
-
-PostOperation {
- { Name Get_Fields ; NameOfPostProcessing PostPro ;
- Operation {
- Print[ p, OnElementsOf VolAll, File "res/p.pos"] ;
- }
- }
-}
+Include "formulations_scalarWaves.pro";
DefineConstant[
C_ = {"-solve -pos -v2", Name "GetDP/9ComputeCommand", Visible 0},
- P_ = {"Get_Fields", Name "GetDP/2PostOperationChoices", Visible 0, ReadOnly 1}
+ R_ = {"WaveEqn_Reso", Name "GetDP/1ResolutionChoices", Visible 1, Choices {"WaveEqn_Reso", "MembraneSys_Reso", "EulerSys_Reso"} },
+ P_ = {"WaveEqn_PostOp", Name "GetDP/2PostOperationChoices", Visible 1, Choices {"WaveEqn_PostOp", "MembraneSys_PostOp", "EulerSys_PostOp"}}
];
-
-
diff --git a/AcademicWaves/scattKite.dat b/AcademicWaves/scattKite.dat
index a4dceb1..c9ff493 100644
--- a/AcademicWaves/scattKite.dat
+++ b/AcademicWaves/scattKite.dat
@@ -3,14 +3,14 @@
// File: parameters
//========================================================
+FLAG_DIM = 2;
+FLAG_RES = RES_FREQ;
//---------------------
// Model choices
//---------------------
DefineConstant[
- Flag_ModelEquations = {1, Choices{1="Wave equation (p)", 2="Wave system (p and u)"},
- Name "Input/01Model equations", Highlight "Black"}
Flag_DomTruncShape = {1, Choices{1="Circular domain", 2="Squared domain"},
Name "Input/1Domain truncation/01Shape of the domain", Highlight "Black"},
Flag_DomTruncMethod_Shape1 = {1, Choices{1="Sommerfeld BC", 2="Bayliss-Turkel BC"},
@@ -25,13 +25,10 @@ If (Flag_DomTruncShape==2)
Flag_DomTruncMethod = Flag_DomTruncMethod_Shape2;
EndIf
-
//---------------------
// Physical parameters
//---------------------
-DIM = 2;
-
// Incident signal
DefineConstant[
@@ -42,7 +39,7 @@ DefineConstant[
theta = { 0., Min 0., Max 2*Pi, Step 0.01,
Name "Input/2Incident signal/12Theta", Highlight "Yellow"}
];
-
+FREQ = k/(2*Pi);
// Geometry -- Fig. 7, page 26 (Ecevit, Reitich, "Multiple scattering iterations for high frequency")
@@ -67,15 +64,7 @@ DefineConstant[
// Numerical parameters
//---------------------
-DefineConstant[
- lc = { 0.05, Min 0.001, Max 1, Step 0.001,
- Name "Input/2Numerical parameters/10Size of a mesh cell", Highlight "Yellow"}
-];
-lc0 = lc;
-lc1 = lc;
-
-ORDER = 1;
-
+LC = lambda/10;
//---------------------
// Geometrical tags
@@ -88,4 +77,3 @@ GAMMA_OBS2 = 202;
DOMAIN = 301;
LAYER = 302;
-
diff --git a/AcademicWaves/scattKite.geo b/AcademicWaves/scattKite.geo
index 4e927fd..02ea46a 100644
--- a/AcademicWaves/scattKite.geo
+++ b/AcademicWaves/scattKite.geo
@@ -3,33 +3,35 @@
// File: GMSH geometry
//========================================================
-Include "scattKite.dat" ;
+Include "scattKite.dat";
+Mesh.CharacteristicLengthMax = LC;
+Mesh.CharacteristicLengthFactor = 1;
+Mesh.Optimize = 1;
//--------------------------------------------------------
// Obstacle 1
//--------------------------------------------------------
-cen0=newp ; Point(newp) = {0,0,0,lc0} ;
-pObs[]+=newp ; Point(newp) = {R0,0,0,lc0} ;
-pObs[]+=newp ; Point(newp) = {0,R1,0,lc0} ;
-pObs[]+=newp ; Point(newp) = {-R0,0,0,lc0} ;
-pObs[]+=newp ; Point(newp) = {0,-R1,0,lc0} ;
-lObs[]+=newl ; Ellipsis(newl) = {pObs[0], cen0, pObs[0], pObs[1]} ;
-lObs[]+=newl ; Ellipsis(newl) = {pObs[1], cen0, pObs[1], pObs[2]} ;
-lObs[]+=newl ; Ellipsis(newl) = {pObs[2], cen0, pObs[2], pObs[3]} ;
-lObs[]+=newl ; Ellipsis(newl) = {pObs[3], cen0, pObs[3], pObs[0]} ;
+cen0=newp; Point(newp) = {0,0,0};
+pObs[]+=newp; Point(newp) = {R0,0,0};
+pObs[]+=newp; Point(newp) = {0,R1,0};
+pObs[]+=newp; Point(newp) = {-R0,0,0};
+pObs[]+=newp; Point(newp) = {0,-R1,0};
+lObs[]+=newl; Ellipsis(newl) = {pObs[0], cen0, pObs[0], pObs[1]};
+lObs[]+=newl; Ellipsis(newl) = {pObs[1], cen0, pObs[1], pObs[2]};
+lObs[]+=newl; Ellipsis(newl) = {pObs[2], cen0, pObs[2], pObs[3]};
+lObs[]+=newl; Ellipsis(newl) = {pObs[3], cen0, pObs[3], pObs[0]};
-llellipsis=newll ; Line Loop(newll) = {lObs[]} ;
-surfellipsis=news ; Plane Surface(news) = {llellipsis} ;
+llellipsis=newll; Line Loop(newll) = {lObs[]};
+surfellipsis=news; Plane Surface(news) = {llellipsis};
-ll[] += llellipsis ;
+ll[] += llellipsis;
Rotate {{0, 0, 1}, {0, 0, 0}, phi} { Surface{surfellipsis}; }
Translate {cenx, ceny, 0} { Surface{surfellipsis}; }
-Physical Line(GAMMA_OBS1) = {-lObs[]} ;
-
+Physical Line(GAMMA_OBS1) = {-lObs[]};
//--------------------------------------------------------
// Obstacle 2 : Kite shaped domain
@@ -37,19 +39,19 @@ Physical Line(GAMMA_OBS1) = {-lObs[]} ;
// .. 0 <= t <= 2*Pi
//--------------------------------------------------------
-t[] = {0:2*Pi:(2*Pi/50)} ;
+
+t[] = {0:2*Pi:(2*Pi/50)};
For kk In {0:#t[]-1}
- pk[]+=newp ; Point(newp) = { factor*(Cos(t[kk])+0.65*Cos(2*t[kk])-0.65), factor*(1.5*Sin(t[kk])), 0., lc0 } ;
+ pk[]+=newp; Point(newp) = { factor*(Cos(t[kk])+0.65*Cos(2*t[kk])-0.65), factor*(1.5*Sin(t[kk])), 0.};
EndFor
-lkite=newl ; Spline(newl) = {pk[],pk[0]} ;
-
-llkite=newll ; Line Loop(newll) = {lkite} ;
-skite=news ; Plane Surface(news) = {newll-1} ;
+lkite=newl; Spline(newl) = {pk[],pk[0]};
-ll[] += llkite ;
+llkite=newll; Line Loop(newll) = {lkite};
+skite=news; Plane Surface(news) = {newll-1};
-Physical Line(GAMMA_OBS2) = {-lkite} ;
+ll[] += llkite;
+Physical Line(GAMMA_OBS2) = {-lkite};
//--------------------------------------------------------
// Main domain and Layers with the boundaries (CIRCLE)
@@ -57,38 +59,38 @@ Physical Line(GAMMA_OBS2) = {-lkite} ;
If (Flag_DomTruncShape==1)
- cen0 = newp ; Point(newp) = {0,0,0,lc0} ;
- pi[]+=newp ; Point(newp) = {Rdom,0,0,lc1} ;
- pi[]+=newp ; Point(newp) = {0,Rdom,0,lc1} ;
- pi[]+=newp ; Point(newp) = {-Rdom,0,0,lc1} ;
- pi[]+=newp ; Point(newp) = {0,-Rdom,0,lc1} ;
- li[]+=newl ; Circle(newl) = {pi[0], cen0, pi[1]} ;
- li[]+=newl ; Circle(newl) = {pi[1], cen0, pi[2]} ;
- li[]+=newl ; Circle(newl) = {pi[2], cen0, pi[3]} ;
- li[]+=newl ; Circle(newl) = {pi[3], cen0, pi[0]} ;
+ cen0 = newp; Point(newp) = {0,0,0};
+ pi[]+=newp; Point(newp) = {Rdom,0,0};
+ pi[]+=newp; Point(newp) = {0,Rdom,0};
+ pi[]+=newp; Point(newp) = {-Rdom,0,0};
+ pi[]+=newp; Point(newp) = {0,-Rdom,0};
+ li[]+=newl; Circle(newl) = {pi[0], cen0, pi[1]};
+ li[]+=newl; Circle(newl) = {pi[1], cen0, pi[2]};
+ li[]+=newl; Circle(newl) = {pi[2], cen0, pi[3]};
+ li[]+=newl; Circle(newl) = {pi[3], cen0, pi[0]};
- llOmega=newll ; Line Loop(llOmega) = li[];
- surfOmega=news ; Plane Surface(news) = {llOmega, ll[]};
+ llOmega=newll; Line Loop(llOmega) = li[];
+ surfOmega=news; Plane Surface(news) = {llOmega, ll[]};
- Physical Line(GAMMA_DOM) = li[] ;
+ Physical Line(GAMMA_DOM) = li[];
Physical Surface(DOMAIN) = {surfOmega};
If (Flag_DomTruncMethod==9)
- cen0 = newp ; Point(newp) = {0,0,0,lc0} ;
- pe[]+=newp ; Point(newp) = {Rdom+Lpml,0,0,lc1} ;
- pe[]+=newp ; Point(newp) = {0,Rdom+Lpml,0,lc1} ;
- pe[]+=newp ; Point(newp) = {-Rdom-Lpml,0,0,lc1} ;
- pe[]+=newp ; Point(newp) = {0,-Rdom-Lpml,0,lc1} ;
- le[]+=newl ; Circle(newl) = {pe[0], cen0, pe[1]} ;
- le[]+=newl ; Circle(newl) = {pe[1], cen0, pe[2]} ;
- le[]+=newl ; Circle(newl) = {pe[2], cen0, pe[3]} ;
- le[]+=newl ; Circle(newl) = {pe[3], cen0, pe[0]} ;
+ cen0 = newp; Point(newp) = {0,0,0};
+ pe[]+=newp; Point(newp) = {Rdom+Lpml,0,0};
+ pe[]+=newp; Point(newp) = {0,Rdom+Lpml,0};
+ pe[]+=newp; Point(newp) = {-Rdom-Lpml,0,0};
+ pe[]+=newp; Point(newp) = {0,-Rdom-Lpml,0};
+ le[]+=newl; Circle(newl) = {pe[0], cen0, pe[1]};
+ le[]+=newl; Circle(newl) = {pe[1], cen0, pe[2]};
+ le[]+=newl; Circle(newl) = {pe[2], cen0, pe[3]};
+ le[]+=newl; Circle(newl) = {pe[3], cen0, pe[0]};
- llLayer=newll ; Line Loop(llLayer) = le[];
- surfLayer=news ; Plane Surface(news) = {llOmega, llLayer};
+ llLayer=newll; Line Loop(llLayer) = le[];
+ surfLayer=news; Plane Surface(news) = {llOmega, llLayer};
- Physical Line(GAMMA_LAY) = le[] ;
+ Physical Line(GAMMA_LAY) = le[];
Physical Surface(LAYER) = {surfLayer};
EndIf
@@ -102,37 +104,37 @@ EndIf
If (Flag_DomTruncShape==2)
- cen0=newp ; Point(newp) = {0,0,0,lc0} ;
- pi[]+=newp ; Point(newp) = { Ldom/2, Ldom/2,0,lc1} ;
- pi[]+=newp ; Point(newp) = {-Ldom/2, Ldom/2,0,lc1} ;
- pi[]+=newp ; Point(newp) = {-Ldom/2,-Ldom/2,0,lc1} ;
- pi[]+=newp ; Point(newp) = { Ldom/2,-Ldom/2,0,lc1} ;
- li[]+=newl ; Line(newl) = {pi[0], pi[1]} ;
- li[]+=newl ; Line(newl) = {pi[1], pi[2]} ;
- li[]+=newl ; Line(newl) = {pi[2], pi[3]} ;
- li[]+=newl ; Line(newl) = {pi[3], pi[0]} ;
+ cen0=newp; Point(newp) = {0,0,0};
+ pi[]+=newp; Point(newp) = { Ldom/2, Ldom/2,0};
+ pi[]+=newp; Point(newp) = {-Ldom/2, Ldom/2,0};
+ pi[]+=newp; Point(newp) = {-Ldom/2,-Ldom/2,0};
+ pi[]+=newp; Point(newp) = { Ldom/2,-Ldom/2,0};
+ li[]+=newl; Line(newl) = {pi[0], pi[1]};
+ li[]+=newl; Line(newl) = {pi[1], pi[2]};
+ li[]+=newl; Line(newl) = {pi[2], pi[3]};
+ li[]+=newl; Line(newl) = {pi[3], pi[0]};
- llOmega=newll ; Line Loop(llOmega) = li[];
- surfOmega=news ; Plane Surface(news) = {llOmega, ll[]};
+ llOmega=newll; Line Loop(llOmega) = li[];
+ surfOmega=news; Plane Surface(news) = {llOmega, ll[]};
- Physical Line(GAMMA_DOM) = li[] ;
+ Physical Line(GAMMA_DOM) = li[];
Physical Surface(DOMAIN) = {surfOmega};
If (Flag_DomTruncMethod==9)
- pe[]+=newp ; Point(newp) = { Ldom/2+Lpml, Ldom/2+Lpml,0,lc1} ;
- pe[]+=newp ; Point(newp) = {-Ldom/2-Lpml, Ldom/2+Lpml,0,lc1} ;
- pe[]+=newp ; Point(newp) = {-Ldom/2-Lpml,-Ldom/2-Lpml,0,lc1} ;
- pe[]+=newp ; Point(newp) = { Ldom/2+Lpml,-Ldom/2-Lpml,0,lc1} ;
- le[]+=newl ; Line(newl) = {pe[0], pe[1]} ;
- le[]+=newl ; Line(newl) = {pe[1], pe[2]} ;
- le[]+=newl ; Line(newl) = {pe[2], pe[3]} ;
- le[]+=newl ; Line(newl) = {pe[3], pe[0]} ;
+ pe[]+=newp; Point(newp) = { Ldom/2+Lpml, Ldom/2+Lpml,0};
+ pe[]+=newp; Point(newp) = {-Ldom/2-Lpml, Ldom/2+Lpml,0};
+ pe[]+=newp; Point(newp) = {-Ldom/2-Lpml,-Ldom/2-Lpml,0};
+ pe[]+=newp; Point(newp) = { Ldom/2+Lpml,-Ldom/2-Lpml,0};
+ le[]+=newl; Line(newl) = {pe[0], pe[1]};
+ le[]+=newl; Line(newl) = {pe[1], pe[2]};
+ le[]+=newl; Line(newl) = {pe[2], pe[3]};
+ le[]+=newl; Line(newl) = {pe[3], pe[0]};
- llLayer=newll ; Line Loop(llLayer) = le[];
- surfLayer=news ; Plane Surface(news) = {llOmega, llLayer};
+ llLayer=newll; Line Loop(llLayer) = le[];
+ surfLayer=news; Plane Surface(news) = {llOmega, llLayer};
- Physical Line(GAMMA_LAY) = le[] ;
+ Physical Line(GAMMA_LAY) = le[];
Physical Surface(LAYER) = {surfLayer};
EndIf
diff --git a/AcademicWaves/scattKite.pro b/AcademicWaves/scattKite.pro
index 2ba7374..796a5ad 100644
--- a/AcademicWaves/scattKite.pro
+++ b/AcademicWaves/scattKite.pro
@@ -3,113 +3,86 @@
// File: GetDP simulation
//========================================================
-Include "scattKite.dat" ;
+Include "scattKite.dat";
Group {
- Omega = Region[{DOMAIN}] ;
- GammaD = Region[{GAMMA_OBS1, GAMMA_OBS2}] ;
+ Omega = Region[{DOMAIN}];
+ GammaD = Region[{GAMMA_OBS1, GAMMA_OBS2}];
If(Flag_DomTruncMethod==1)
- GammaR = Region[{GAMMA_DOM}] ;
+ GammaR = Region[{GAMMA_DOM}];
EndIf
If(Flag_DomTruncMethod==2)
- GammaBT = Region[{GAMMA_DOM}] ;
+ GammaBT = Region[{GAMMA_DOM}];
EndIf
If(Flag_DomTruncMethod==9)
- OmegaPml = Region[{LAYER}] ;
- GammaD += Region[{GAMMA_LAY}] ;
+ OmegaPml = Region[{LAYER}];
+ GammaD += Region[{GAMMA_LAY}];
EndIf
}
Function {
-
- I[] = Complex[0, 1] ;
-
- Freq = 1./lambda;
- kVec[] = Vector[ k*Cos[theta], k*Sin[theta], 0] ;
- pInc[] = Complex[ Cos[kVec[]*XYZ[]], Sin[kVec[]*XYZ[]] ] ;
- uInc[] = -pInc[] * kVec[]/k ;
-
+ I[] = Complex[0, 1];
+
+ // Incident signal
+ kVec[] = Vector[ k*Cos[theta], k*Sin[theta], 0];
+ pInc[] = Complex[ Cos[kVec[]*XYZ[]], Sin[kVec[]*XYZ[]] ];
+
// BC: Bayliss-Turkel
If(Flag_DomTruncMethod==2)
- alphaBT[] = 1/(2*Ldom) - I[]/(8*k*Ldom^2*(1+I[]/(k*Ldom))) ;
- betaBT[] = - 1/(2*I[]*k*(1+I[]/(k*Ldom))) ;
+ alphaBT[] = 1/(2*Ldom) - I[]/(8*k*Ldom^2*(1+I[]/(k*Ldom)));
+ betaBT[] = - 1/(2*I[]*k*(1+I[]/(k*Ldom)));
EndIf
// PML
If(Flag_DomTruncMethod==9)
xLoc[] = Fabs[X[]]-Ldom/2;
yLoc[] = Fabs[Y[]]-Ldom/2;
- absFuncX[] = (xLoc[]>=0) ? 1/(Lpml-xLoc[]) : 0 ;
- absFuncY[] = (yLoc[]>=0) ? 1/(Lpml-yLoc[]) : 0 ;
+ absFuncX[] = (xLoc[]>=0) ? Log[Lpml-xLoc[]] : 0;
+ absFuncY[] = (yLoc[]>=0) ? Log[Lpml-yLoc[]] : 0;
hx[] = Complex[1, absFuncX[]/k];
hy[] = Complex[1, absFuncY[]/k];
- layerScal[OmegaPml] = 1/(hx[]*hy[]);
- layerTens[OmegaPml] = TensorDiag[hy[]/hx[], hx[]/hy[], 1.];
+ pmlScal[OmegaPml] = 1/(hx[]*hy[]);
+ pmlTens[OmegaPml] = TensorDiag[hy[]/hx[], hx[]/hy[], 1.];
EndIf
-
}
Constraint {
- { Name pConstraint ;
+ { Name pConstraint;
Case {
- { Region Region[{GAMMA_OBS1, GAMMA_OBS2}] ; Value -pInc[] ; }
+ { Region Region[{GAMMA_OBS1, GAMMA_OBS2}]; Value -pInc[]; }
If(Flag_DomTruncMethod==9)
- { Region Region[{GAMMA_LAY}] ; Value 0. ; }
+ { Region Region[{GAMMA_LAY}]; Value 0.; }
EndIf
}
}
- { Name uConstraint ;
- Case {
- { }
- }
- }
+ { Name uConstraint ; Case {} }
}
-Include "formulations_scalarWaves.pro" ;
+Include "formulations_scalarWaves.pro";
PostProcessing {
- If (Flag_ModelEquations==1)
- { Name PostPro ; NameOfFormulation WaveEquation_Form ;
+ { Name WaveEqn2_PostPro; NameOfFormulation WaveEqn_Form;
Quantity {
- { Name pSct ; Value{ Local{ [ {p} ] ; In VolAll ; Jacobian JVol ; } } }
- { Name pInc ; Value{ Local{ [ pInc[] ] ; In VolAll ; Jacobian JVol ; } } }
- { Name pTot ; Value{ Local{ [ pInc[]+{p} ] ; In VolAll ; Jacobian JVol ; } } }
+ { Name pSct; Value{ Local{ [ {p} ]; In VolAll; Jacobian JVol; }}}
+ { Name pInc; Value{ Local{ [ pInc[] ]; In VolAll; Jacobian JVol; }}}
+ { Name pTot; Value{ Local{ [ pInc[]+{p} ]; In VolAll; Jacobian JVol; }}}
}
}
- EndIf
- If (Flag_ModelEquations==2)
- { Name PostPro ; NameOfFormulation EulerEquations_Form ;
- Quantity {
- { Name pSct ; Value{ Local{ [ {p} ] ; In VolAll ; Jacobian JVol ; } } }
- { Name pInc ; Value{ Local{ [ pInc[] ] ; In VolAll ; Jacobian JVol ; } } }
- { Name pTot ; Value{ Local{ [ pInc[]+{p} ] ; In VolAll ; Jacobian JVol ; } } }
- { Name uSct ; Value{ Local{ [ {u} ] ; In VolAll ; Jacobian JVol ; } } }
- { Name uInc ; Value{ Local{ [ uInc[] ] ; In VolAll ; Jacobian JVol ; } } }
- { Name uTot ; Value{ Local{ [ uInc[]+{u} ] ; In VolAll ; Jacobian JVol ; } } }
- }
- }
- EndIf
}
PostOperation {
- { Name Get_Fields ; NameOfPostProcessing PostPro ;
+ { Name WaveEqn2_PostOp; NameOfPostProcessing WaveEqn2_PostPro;
Operation {
- Print[ pSct, OnElementsOf VolAll, File "res/pSct.pos"] ;
- //Print[ pInc, OnElementsOf VolAll, File "res/pInc.pos"] ;
- //Print[ pTot, OnElementsOf VolAll, File "res/pTot.pos"] ;
- If (Flag_ModelEquations==2)
- //Print[ uSct, OnElementsOf VolAll, File "res/uSct.pos"] ;
- //Print[ uInc, OnElementsOf VolAll, File "res/uInc.pos"] ;
- //Print[ uTot, OnElementsOf VolAll, File "res/uTot.pos"] ;
- EndIf
+ Print[ pSct, OnElementsOf VolAll, File "res/pSct.pos"];
+ Print[ pInc, OnElementsOf VolAll, File "res/pInc.pos"];
+ Print[ pTot, OnElementsOf VolAll, File "res/pTot.pos"];
}
}
}
DefineConstant[
C_ = {"-solve -pos -v2", Name "GetDP/9ComputeCommand", Visible 0},
- P_ = {"Get_Fields", Name "GetDP/2PostOperationChoices", Visible 0, ReadOnly 1}
+ R_ = {"WaveEqn_Reso", Name "GetDP/1ResolutionChoices", Visible 0, Choices {"WaveEqn_Reso"} },
+ P_ = {"WaveEqn2_PostOp", Name "GetDP/2PostOperationChoices", Visible 0, Choices {"WaveEqn2_PostOp"}}
];
-
-
diff --git a/AcademicWaves/scattPlates.dat b/AcademicWaves/scattPlates.dat
index b050733..53fa31a 100644
--- a/AcademicWaves/scattPlates.dat
+++ b/AcademicWaves/scattPlates.dat
@@ -3,14 +3,14 @@
// File: parameters
//========================================================
+FLAG_DIM = 2;
+FLAG_RES = RES_FREQ;
//---------------------
// Model choices
//---------------------
DefineConstant[
- Flag_ModelEquations = {1, Choices{1="Wave equation (p)", 2="Wave system (p and u)"},
- Name "Input/01Model equations", Highlight "Black"}
Flag_DomTruncShape = {1, Choices{1="Circular domain", 2="Squared domain"},
Name "Input/1Domain truncation/01Shape of the domain", Highlight "Black"},
Flag_DomTruncMethod_Shape1 = {1, Choices{1="Sommerfeld BC", 2="Bayliss-Turkel BC"},
@@ -40,6 +40,7 @@ DefineConstant[
theta = { Pi/4, Min 0., Max 2*Pi, Step 0.01,
Name "Input/2Incident signal/12Theta", Highlight "Yellow"}
];
+FREQ = k/(2*Pi);
// Geometry
@@ -57,20 +58,11 @@ DefineConstant[
Name "Input/1Domain truncation/12Thickness of the layer", Highlight "Yellow", Visible (Flag_DomTruncMethod==9)}
];
-
//---------------------
// Numerical parameters
//---------------------
-DefineConstant[
- lc = { 0.01, Min 0.001, Max 1, Step 0.001,
- Name "Input/2Numerical parameters/10Size of a mesh cell", Highlight "Yellow"}
-];
-lc0 = lc;
-lc1 = lc;
-
-Flag_SecondOrder = 0;
-
+LC = lambda/10;
//---------------------
// Geometrical tags
diff --git a/AcademicWaves/scattPlates.geo b/AcademicWaves/scattPlates.geo
index 34e7091..56faf98 100644
--- a/AcademicWaves/scattPlates.geo
+++ b/AcademicWaves/scattPlates.geo
@@ -5,15 +5,18 @@
Include "scattPlates.dat";
+Mesh.CharacteristicLengthMax = LC;
+Mesh.CharacteristicLengthFactor = 1;
+Mesh.Optimize = 1;
//--------------------------------------------------------
// Obstacle 1
//--------------------------------------------------------
-pObs1[]+=newp ; Point(newp) = { dx/2,-dy/2+dist, 0, lc0} ;
-pObs1[]+=newp ; Point(newp) = { dx/2, dy/2+dist, 0, lc0} ;
-pObs1[]+=newp ; Point(newp) = {-dx/2, dy/2+dist, 0, lc0} ;
-pObs1[]+=newp ; Point(newp) = {-dx/2,-dy/2+dist, 0, lc0} ;
+pObs1[]+=newp ; Point(newp) = { dx/2,-dy/2+dist, 0} ;
+pObs1[]+=newp ; Point(newp) = { dx/2, dy/2+dist, 0} ;
+pObs1[]+=newp ; Point(newp) = {-dx/2, dy/2+dist, 0} ;
+pObs1[]+=newp ; Point(newp) = {-dx/2,-dy/2+dist, 0} ;
lObs1[]+=newl ; Line(newl) = {pObs1[0],pObs1[1]} ;
lObs1[]+=newl ; Line(newl) = {pObs1[1],pObs1[2]} ;
@@ -30,10 +33,10 @@ Physical Line(GAMMA_OBS1) = {-lObs1[]} ;
// Obstacle 2
//--------------------------------------------------------
-pObs2[]+=newp ; Point(newp) = { dx/2,-dy/2-dist, 0, lc0} ;
-pObs2[]+=newp ; Point(newp) = { dx/2, dy/2-dist, 0, lc0} ;
-pObs2[]+=newp ; Point(newp) = {-dx/2, dy/2-dist, 0, lc0} ;
-pObs2[]+=newp ; Point(newp) = {-dx/2,-dy/2-dist, 0, lc0} ;
+pObs2[]+=newp ; Point(newp) = { dx/2,-dy/2-dist, 0} ;
+pObs2[]+=newp ; Point(newp) = { dx/2, dy/2-dist, 0} ;
+pObs2[]+=newp ; Point(newp) = {-dx/2, dy/2-dist, 0} ;
+pObs2[]+=newp ; Point(newp) = {-dx/2,-dy/2-dist, 0} ;
lObs2[]+=newl ; Line(newl) = {pObs2[0],pObs2[1]} ;
lObs2[]+=newl ; Line(newl) = {pObs2[1],pObs2[2]} ;
@@ -52,11 +55,11 @@ Physical Line(GAMMA_OBS2) = {-lObs2[]} ;
If (Flag_DomTruncShape==1)
- cen0 = newp ; Point(newp) = {0,0,0,lc0} ;
- pi[]+=newp ; Point(newp) = {Rdom,0,0,lc1} ;
- pi[]+=newp ; Point(newp) = {0,Rdom,0,lc1} ;
- pi[]+=newp ; Point(newp) = {-Rdom,0,0,lc1} ;
- pi[]+=newp ; Point(newp) = {0,-Rdom,0,lc1} ;
+ cen0 = newp ; Point(newp) = {0,0,0} ;
+ pi[]+=newp ; Point(newp) = {Rdom,0,0} ;
+ pi[]+=newp ; Point(newp) = {0,Rdom,0} ;
+ pi[]+=newp ; Point(newp) = {-Rdom,0,0} ;
+ pi[]+=newp ; Point(newp) = {0,-Rdom,0} ;
li[]+=newl ; Circle(newl) = {pi[0], cen0, pi[1]} ;
li[]+=newl ; Circle(newl) = {pi[1], cen0, pi[2]} ;
li[]+=newl ; Circle(newl) = {pi[2], cen0, pi[3]} ;
@@ -70,11 +73,11 @@ If (Flag_DomTruncShape==1)
If (Flag_DomTruncMethod==9)
- cen0 = newp ; Point(newp) = {0,0,0,lc0} ;
- pe[]+=newp ; Point(newp) = {Rdom+Lpml,0,0,lc1} ;
- pe[]+=newp ; Point(newp) = {0,Rdom+Lpml,0,lc1} ;
- pe[]+=newp ; Point(newp) = {-Rdom-Lpml,0,0,lc1} ;
- pe[]+=newp ; Point(newp) = {0,-Rdom-Lpml,0,lc1} ;
+ cen0 = newp ; Point(newp) = {0,0,0} ;
+ pe[]+=newp ; Point(newp) = {Rdom+Lpml,0,0} ;
+ pe[]+=newp ; Point(newp) = {0,Rdom+Lpml,0} ;
+ pe[]+=newp ; Point(newp) = {-Rdom-Lpml,0,0} ;
+ pe[]+=newp ; Point(newp) = {0,-Rdom-Lpml,0} ;
le[]+=newl ; Circle(newl) = {pe[0], cen0, pe[1]} ;
le[]+=newl ; Circle(newl) = {pe[1], cen0, pe[2]} ;
le[]+=newl ; Circle(newl) = {pe[2], cen0, pe[3]} ;
@@ -97,11 +100,11 @@ EndIf
If (Flag_DomTruncShape==2)
- cen0=newp ; Point(newp) = {0,0,0,lc0} ;
- pi[]+=newp ; Point(newp) = { Ldom/2, Ldom/2,0,lc1} ;
- pi[]+=newp ; Point(newp) = {-Ldom/2, Ldom/2,0,lc1} ;
- pi[]+=newp ; Point(newp) = {-Ldom/2,-Ldom/2,0,lc1} ;
- pi[]+=newp ; Point(newp) = { Ldom/2,-Ldom/2,0,lc1} ;
+ cen0=newp ; Point(newp) = {0,0,0} ;
+ pi[]+=newp ; Point(newp) = { Ldom/2, Ldom/2,0} ;
+ pi[]+=newp ; Point(newp) = {-Ldom/2, Ldom/2,0} ;
+ pi[]+=newp ; Point(newp) = {-Ldom/2,-Ldom/2,0} ;
+ pi[]+=newp ; Point(newp) = { Ldom/2,-Ldom/2,0} ;
li[]+=newl ; Line(newl) = {pi[0], pi[1]} ;
li[]+=newl ; Line(newl) = {pi[1], pi[2]} ;
li[]+=newl ; Line(newl) = {pi[2], pi[3]} ;
@@ -115,10 +118,10 @@ If (Flag_DomTruncShape==2)
If (Flag_DomTruncMethod==9)
- pe[]+=newp ; Point(newp) = { Ldom/2+Lpml, Ldom/2+Lpml,0,lc1} ;
- pe[]+=newp ; Point(newp) = {-Ldom/2-Lpml, Ldom/2+Lpml,0,lc1} ;
- pe[]+=newp ; Point(newp) = {-Ldom/2-Lpml,-Ldom/2-Lpml,0,lc1} ;
- pe[]+=newp ; Point(newp) = { Ldom/2+Lpml,-Ldom/2-Lpml,0,lc1} ;
+ pe[]+=newp ; Point(newp) = { Ldom/2+Lpml, Ldom/2+Lpml,0} ;
+ pe[]+=newp ; Point(newp) = {-Ldom/2-Lpml, Ldom/2+Lpml,0} ;
+ pe[]+=newp ; Point(newp) = {-Ldom/2-Lpml,-Ldom/2-Lpml,0} ;
+ pe[]+=newp ; Point(newp) = { Ldom/2+Lpml,-Ldom/2-Lpml,0} ;
le[]+=newl ; Line(newl) = {pe[0], pe[1]} ;
le[]+=newl ; Line(newl) = {pe[1], pe[2]} ;
le[]+=newl ; Line(newl) = {pe[2], pe[3]} ;
diff --git a/AcademicWaves/scattPlates.pro b/AcademicWaves/scattPlates.pro
index f2eebcc..80a1bfc 100644
--- a/AcademicWaves/scattPlates.pro
+++ b/AcademicWaves/scattPlates.pro
@@ -3,108 +3,86 @@
// File: GetDP simulation
//========================================================
-Include "scattPlates.dat" ;
+Include "scattPlates.dat";
Group {
- Omega = Region[{DOMAIN}] ;
- GammaD = Region[{GAMMA_OBS1, GAMMA_OBS2}] ;
-
+ Omega = Region[{DOMAIN}];
+ GammaD = Region[{GAMMA_OBS1, GAMMA_OBS2}];
+
If(Flag_DomTruncMethod==1)
- GammaR = Region[{GAMMA_DOM}] ;
+ GammaR = Region[{GAMMA_DOM}];
EndIf
If(Flag_DomTruncMethod==2)
- GammaBT = Region[{GAMMA_DOM}] ;
+ GammaBT = Region[{GAMMA_DOM}];
EndIf
If(Flag_DomTruncMethod==9)
- OmegaPml = Region[{LAYER}] ;
- GammaD += Region[{GAMMA_LAY}] ;
+ OmegaPml = Region[{LAYER}];
+ GammaD += Region[{GAMMA_LAY}];
EndIf
}
Function {
-
- I[] = Complex[0, 1] ;
-
+ I[] = Complex[0, 1];
+
// Incident signal
- kVec[] = Vector[ k*Cos[theta], k*Sin[theta], 0] ;
- pInc[] = Complex[ Cos[kVec[]*XYZ[]], Sin[kVec[]*XYZ[]] ] ;
-
+ kVec[] = Vector[ k*Cos[theta], k*Sin[theta], 0];
+ pInc[] = Complex[ Cos[kVec[]*XYZ[]], Sin[kVec[]*XYZ[]] ];
+
// BC: Bayliss-Turkel
If(Flag_DomTruncMethod==2)
- alphaBT[] = 1/(2*Ldom) - I[]/(8*k*Ldom^2*(1+I[]/(k*Ldom))) ;
- betaBT[] = - 1/(2*I[]*k*(1+I[]/(k*Ldom))) ;
+ alphaBT[] = 1/(2*Ldom) - I[]/(8*k*Ldom^2*(1+I[]/(k*Ldom)));
+ betaBT[] = - 1/(2*I[]*k*(1+I[]/(k*Ldom)));
EndIf
-
+
// PML
If(Flag_DomTruncMethod==9)
xLoc[] = Fabs[X[]]-Ldom/2;
yLoc[] = Fabs[Y[]]-Ldom/2;
- absFuncX[] = (xLoc[]>=0) ? 1/(Lpml-xLoc[]) : 0 ;
- absFuncY[] = (yLoc[]>=0) ? 1/(Lpml-yLoc[]) : 0 ;
+ absFuncX[] = (xLoc[]>=0) ? Log[Lpml-xLoc[]] : 0;
+ absFuncY[] = (yLoc[]>=0) ? Log[Lpml-yLoc[]] : 0;
hx[] = Complex[1, absFuncX[]/k];
hy[] = Complex[1, absFuncY[]/k];
- layerScal[OmegaPml] = 1/(hx[]*hy[]);
- layerTens[OmegaPml] = TensorDiag[hy[]/hx[], hx[]/hy[], 1.];
+ pmlScal[OmegaPml] = 1/(hx[]*hy[]);
+ pmlTens[OmegaPml] = TensorDiag[hy[]/hx[], hx[]/hy[], 1.];
EndIf
-
}
Constraint {
- { Name pConstraint ;
+ { Name pConstraint;
Case {
- { Region Region[{GAMMA_OBS1, GAMMA_OBS2}] ; Value -pInc[] ; }
+ { Region Region[{GAMMA_OBS1, GAMMA_OBS2}]; Value -pInc[]; }
If(Flag_DomTruncMethod==9)
- { Region Region[{GAMMA_LAY}] ; Value 0. ; }
+ { Region Region[{GAMMA_LAY}]; Value 0.; }
EndIf
}
}
- { Name uConstraint ;
- Case {
- }
- }
+ { Name uConstraint ; Case {} }
}
-If (Flag_ModelEquations==1)
- Include "form_scalarWaveEqn_frequency.pro" ;
-EndIf
-If (Flag_ModelEquations==2)
- Include "form_scalarWaveSyst2D_frequency.pro" ;
-EndIf
+Include "formulations_scalarWaves.pro";
PostProcessing {
- { Name PostPro ; NameOfFormulation Form ;
+ { Name WaveEqn2_PostPro; NameOfFormulation WaveEqn_Form;
Quantity {
- { Name pSct ; Value{ Local{ [ {p} ] ; In VolAll ; Jacobian JVol ; } } }
- { Name pInc ; Value{ Local{ [ pInc[] ] ; In VolAll ; Jacobian JVol ; } } }
- { Name pTot ; Value{ Local{ [ pInc[]+{p} ] ; In VolAll ; Jacobian JVol ; } } }
- If (Flag_ModelEquations==2)
- { Name uSct ; Value{ Local{ [ {u} ] ; In VolAll ; Jacobian JVol ; } } }
- { Name uInc ; Value{ Local{ [ uInc[] ] ; In VolAll ; Jacobian JVol ; } } }
- { Name uTot ; Value{ Local{ [ uInc[]+{u} ] ; In VolAll ; Jacobian JVol ; } } }
- EndIf
+ { Name pSct; Value{ Local{ [ {p} ]; In VolAll; Jacobian JVol; }}}
+ { Name pInc; Value{ Local{ [ pInc[] ]; In VolAll; Jacobian JVol; }}}
+ { Name pTot; Value{ Local{ [ pInc[]+{p} ]; In VolAll; Jacobian JVol; }}}
}
}
}
PostOperation {
- { Name Get_Fields ; NameOfPostProcessing PostPro ;
+ { Name WaveEqn2_PostOp; NameOfPostProcessing WaveEqn2_PostPro;
Operation {
- Print[ pSct, OnElementsOf VolAll, File "res/pSct.pos"] ;
- Print[ pInc, OnElementsOf VolAll, File "res/pInc.pos"] ;
- Print[ pTot, OnElementsOf VolAll, File "res/pTot.pos"] ;
- If (Flag_ModelEquations==2)
- Print[ uSct, OnElementsOf VolAll, File "res/uSct.pos"] ;
- Print[ uInc, OnElementsOf VolAll, File "res/uInc.pos"] ;
- Print[ uTot, OnElementsOf VolAll, File "res/uTot.pos"] ;
- EndIf
+ Print[ pSct, OnElementsOf VolAll, File "res/pSct.pos"];
+ Print[ pInc, OnElementsOf VolAll, File "res/pInc.pos"];
+ Print[ pTot, OnElementsOf VolAll, File "res/pTot.pos"];
}
}
}
DefineConstant[
- R_ = {"Reso", Name "GetDP/1ResolutionChoices", Visible 0},
C_ = {"-solve -pos -v2", Name "GetDP/9ComputeCommand", Visible 0},
- P_ = {"Get_Fields", Name "GetDP/2PostOperationChoices", Visible 0, ReadOnly 1}
+ R_ = {"WaveEqn_Reso", Name "GetDP/1ResolutionChoices", Visible 0, Choices {"WaveEqn_Reso"} },
+ P_ = {"WaveEqn2_PostOp", Name "GetDP/2PostOperationChoices", Visible 0, Choices {"WaveEqn2_PostOp"}}
];
-
-
diff --git a/AcademicWaves/scattSphere.dat b/AcademicWaves/scattSphere.dat
index 62350e3..cb14909 100644
--- a/AcademicWaves/scattSphere.dat
+++ b/AcademicWaves/scattSphere.dat
@@ -3,6 +3,8 @@
// File: parameters
//========================================================
+FLAG_DIM = 3;
+FLAG_RES = RES_FREQ;
//---------------------
// Physical parameters
@@ -26,8 +28,8 @@ z0 = Sqrt[mu0/ep0];
// Incident field
-Freq = 1e8;
-k = 2*Pi*Freq/c0;
+FREQ = 1e8;
+k = 2*Pi*FREQ/c0;
// Domain truncation
@@ -35,13 +37,11 @@ Flag_DomTruncMethod = 9;
// 1: Sommerfeld
// 9: PML
-
//---------------------
// Numerical parameters
//---------------------
-lc = 0.3;
-
+LC = 0.3;
//---------------------
// Geometrical tags
diff --git a/AcademicWaves/scattSphere.geo b/AcademicWaves/scattSphere.geo
index a8c0e47..c190e1b 100644
--- a/AcademicWaves/scattSphere.geo
+++ b/AcademicWaves/scattSphere.geo
@@ -5,13 +5,12 @@
Include "scattSphere.dat";
-Mesh.CharacteristicLengthMax = lc;
+Mesh.CharacteristicLengthMax = LC;
Mesh.CharacteristicLengthFactor = 1;
Mesh.Optimize = 1;
c0 = newp ; Point(newp) = {0,0,0};
-
// SPHERE
ps0[]+=newp ; Point(newp) = { R0, 0, 0};
diff --git a/AcademicWaves/scattSphere.pro b/AcademicWaves/scattSphere.pro
index 058ba33..ffe85be 100644
--- a/AcademicWaves/scattSphere.pro
+++ b/AcademicWaves/scattSphere.pro
@@ -63,9 +63,9 @@ Function {
DampingProfileY[] = (Fabs[Y[]]>=ay) ? c0 * 1/(lpml-(Fabs[Y[]]-ay)) : 0 ;
DampingProfileZ[] = (Fabs[Z[]]>=az) ? c0 * 1/(lpml-(Fabs[Z[]]-az)) : 0 ;
- cX[] = Complex[1, -1/(2*Pi*Freq) * DampingProfileX[]] ;
- cY[] = Complex[1, -1/(2*Pi*Freq) * DampingProfileY[]] ;
- cZ[] = Complex[1, -1/(2*Pi*Freq) * DampingProfileZ[]] ;
+ cX[] = Complex[1, -1/(2*Pi*FREQ) * DampingProfileX[]] ;
+ cY[] = Complex[1, -1/(2*Pi*FREQ) * DampingProfileY[]] ;
+ cZ[] = Complex[1, -1/(2*Pi*FREQ) * DampingProfileZ[]] ;
epsilon[Omega_Pml] = ep0 * TensorDiag[(cY[]*cZ[])/cX[], (cX[]*cZ[])/cY[], (cX[]*cY[])/cZ[]] ;
nu[Omega_Pml] = nu0 * TensorDiag[cX[]/(cY[]*cZ[]), cY[]/(cX[]*cZ[]), cZ[]/(cX[]*cY[])] ;
@@ -76,7 +76,7 @@ Function {
}
-Include "form_vectorWaveEqn_frequency.pro" ;
+Include "formulations_vectorWaves.pro" ;
PostProcessing {
{ Name PostPro ; NameOfFormulation Form ;
--
GitLab