maxwell.pro

Group{
  Omega = Region[1]; // Domain
  Bndry = Region[2]; // Boundary
}

Function{
  // Physical data //
  DefineConstant[angularFreqRe = 5,  // [rad/s]
                 angularFreqIm = 0]; // [rad/s]
  C0    = 299792458;                 // [m/s]
  aFC[] = Complex[angularFreqRe,
                  angularFreqIm];    // [rad/m]
  k[]   = aFC[]/C0;                  // [rad/m]

  Printf["Angular frequency set to: %g + i*%g [rad/s]",
         angularFreqRe, angularFreqIm];

  // Algebraic data //
  DefineConstant[nRHS = 1];       // Number of RHS for this run
  For i In {0:nRHS-1}
    DefineConstant[x~{i}() = {},  // Solution
                   b~{i}() = {}]; // Right hand side
  EndFor

  // Control data //
  DefineConstant[doInit    = 0,          // Should I initialize some stuff?
                 doSolve   = 0,          // Should I solve Ax = b?
                 doPostpro = 0,          // Should I only create a view for x()?
                 doApply   = 0,          // Should I only apply x(): x <- Ax?
                 fileName  = "eig.pos"]; // Postpro file name
}

Jacobian{
  { Name JVol;
    Case{
      { Region All; Jacobian Vol; }
    }
  }
}

Integration{
  { Name IP2;
    Case{
      { Type Gauss;
        Case{
          { GeoElement Line;        NumberOfPoints 2; }
          { GeoElement Triangle;    NumberOfPoints 3; }
          { GeoElement Tetrahedron; NumberOfPoints 4; }
        }
      }
    }
  }
}

Constraint{
  { Name Dirichlet; Type Assign;
    Case{
      { Region Bndry; Value 0; }
    }
  }
}

FunctionSpace{
  { Name HCurl; Type Form1;
    BasisFunction{ // Nedelec
      { Name se; NameOfCoef ee; Function BF_Edge;
        Support Omega; Entity EdgesOf[All]; }
    }

    Constraint{
      { NameOfCoef ee; EntityType EdgesOf; NameOfConstraint Dirichlet; }
    }
  }
}

Formulation{
  { Name Maxwell; Type FemEquation;
    Quantity{
      { Name e; Type Local; NameOfSpace HCurl; }
    }
    Equation{
      Galerkin{ [         Dof{d e} , {d e}];
        In Omega; Integration IP2; Jacobian JVol; }
      Galerkin{ [-k[]^2 * Dof{  e} , {  e}];
        In Omega; Integration IP2; Jacobian JVol; }
    }
  }
}

Resolution{
  { Name Maxwell;
    System{
      { Name A; NameOfFormulation Maxwell; Type Complex; }
    }
    Operation{
      Generate[A];

      If(doInit)
        CopySolution[A, x~{0}()];
      EndIf

      If(doSolve)
        // Full solve for first RHS
        CopyRightHandSide[b~{0}(), A];
        Solve[A];
        CopySolution[A, x~{0}()];

        // SolveAgain for other RHSs
        For i In {1:nRHS-1}
          CopyRightHandSide[b~{i}(), A];
          SolveAgain[A];
          CopySolution[A, x~{i}()];
        EndFor
      EndIf

      If(doApply)
        CopySolution[x~{0}(), A];
        Apply[A];
        CopySolution[A, x~{0}()];
      EndIf

      If(doPostpro)
        CopySolution[x~{0}(), A];
        PostOperation[Maxwell];
      EndIf
    }
  }
}

PostProcessing{
  { Name Maxwell; NameOfFormulation Maxwell;
    Quantity{
      { Name e; Value{ Local{ [{e}]; In Omega; Jacobian JVol; } } }
    }
  }
}

PostOperation{
  { Name Maxwell; NameOfPostProcessing Maxwell;
    Operation{
      Print[e, OnElementsOf Omega, File fileName];
    }
  }
}