[Getdp] Multiple right hand side for the same stiffness matrix

[Benjamin Vial]

Dear all,


First I would like to thank all getdp and gmsh authors/developers for such amazing pieces of software

and to make them available open source!


I am currently working on topology optimization problems in electromagnetism.

I'll try to keep it simple: I want to maximize square the norm of the electric field in some sub-domain Omega_target,

by changing the permittivity value in a design sub-domain Omega_design.

I defined a formulation "helmoltz_scalar" for a 2D time-harmonic problem (xy plane), TE polarization

for the z component of the E-field, which I call u.

In matrix form :

A*u = b   (1)

where b is a known source term.


Now gradient-based optimization algorithm needs information about the sensitivity of the

objective function g with respect to the design variable epsilon. The most efficient way to obtain it

is by solving an adjoint problem, namely

A*uadj = dg/du     (2)

In my case g=|u|^2 so the second member of the adjoint equation (2) is dg/du = 2*Re[u]


Since the matrix A is the same for both problems, we do not need to assemble it twice

Is there a way to save time and generate only the right hand side for the adjoint equation?

It looks like the function GenerateRHSGroup might be useful (but it is not documented yet), together with SolveAgain.


This somehow similar to this old thread but maybe there is a more straightforward way to do it?

http://onelab.info/pipermail/getdp/2011/001383.html



Here are the Formulation and Resolution parts of my .pro file:


#############################################################################

Formulation{
    {Name helmoltz_scalar; Type FemEquation;
        Quantity {
        { Name u; Type Local; NameOfSpace Hgradu;}
        }
        Equation {
        Galerkin { [k0^2*CompZZ[epsilonr[]]*Dof{u} , {u}];
        In Omega; Jacobian JVol; Integration Int_1;  }
        Galerkin { [-1/TensorDiag[CompYY[mur[]],CompXX[mur[]],CompXX[mur[]]]*Dof{Grad u} , {Grad u}];
        In Omega; Jacobian JVol; Integration Int_1; }
        Galerkin { [source[] , {u}];
        In Omega; Jacobian JVol; Integration Int_1;  }
        }
    }
    {Name adjoint_eq; Type FemEquation;
        Quantity {
        { Name u; Type Local; NameOfSpace Hgradu;}
        { Name uadj; Type Local; NameOfSpace Hgraduadj;}
        }
        Equation {
        Galerkin { [  k0^2*CompZZ[epsilonr[]]*Dof{uadj} , {uadj}];
        In Omega; Jacobian JVol; Integration Int_1;  }
        Galerkin { [-1/TensorDiag[CompYY[mur[]],CompXX[mur[]],CompXX[mur[]]]*Dof{Grad uadj} , {Grad uadj}];
        In Omega; Jacobian JVol; Integration Int_1; }
        Galerkin { [2*Re[{u}], {uadj}];  //
        In Omega_target; Jacobian JVol; Integration Int_1;  }
        }
    }

}


Resolution {
    { Name all
        System {
        { Name S; NameOfFormulation helmoltz_scalar; Type ComplexValue; Frequency Freq;}
        { Name Sadj; NameOfFormulation adjoint_eq; Type ComplexValue; Frequency Freq;}
        }
        Operation {
        Generate[S] ;Solve[S] ;SaveSolution[S] ;
        Generate[Sadj] ;Solve[Sadj] ; SaveSolution[Sadj] ;
        }

    }
}


#############################################################################




Thanks in advance!


All best


Benjamin



---------------------------------------------------------------------------------------------------
Benjamin Vial
Postdoctoral Research Assistant
Antennas and Electromagnetics research group
School of Electronics Engineering and Computer Science
Queen Mary, University of London
http://benjaminvial.org/
http://antennas.eecs.qmul.ac.uk/
http://www.quest-spatial-transformation.org/team/researchers/dr-benjamin-vial/

-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://onelab.info/pipermail/getdp/attachments/20160908/4dbe55f1/attachment.html>