Electric_Field_Methods::Implicit_E_Field Class Reference

#include <implicitE.h>

Inheritance diagram for Electric_Field_Methods::Implicit_E_Field:

Collaboration diagram for Electric_Field_Methods::Implicit_E_Field:

Public Member Functions
	Implicit_E_Field (EMF1D &emf, const double &deltat, const double &deltax)
void	advance (Algorithms::RK3< State1D > *rk, State1D &Y, collisions &coll, VlasovFunctor1D_implicitE_p2 *rkF)
Public Member Functions inherited from Electric_Field_Methods::Efield_Method
virtual	~Efield_Method ()=0

Private Member Functions
void	Ampere (EMF1D &emf)
void	FindDE (EMF1D &emf)

Private Attributes
Current_xyz	JN
Current_xyz	J0
Current_xyz	J_Ex
Current_xyz	J_Ey
Current_xyz	J_Ez
Efield_xyz	EN
Efield_xyz	E0
Efield_xyz	DE
double	dt
complex< double >	idx
int	Nbc
int	szx

Detailed Description

Definition at line 145 of file implicitE.h.

Constructor & Destructor Documentation

Implicit_E_Field()

Electric_Field_Methods::Implicit_E_Field::Implicit_E_Field	(	EMF1D &	emf,
		const double &	deltat,
		const double &	deltax
	)

Definition at line 249 of file implicitE.cpp.

References advance(), Input::Input_List::BoundaryCells, Input::List(), Nbc, Input::Input_List::NxLocal, and szx.

Referenced by Electric_Field_Methods::Efield_Method::~Efield_Method().

                                                                            ://, int tout_start)
//--------------------------------------------------------------
//  Constructor for the implicit electric field
//--------------------------------------------------------------
         JN(emf),   J0(emf),                     // New current, Current from E = 0
         J_Ex(emf), J_Ey(emf), J_Ez(emf),        // Current due to the effect of Ex, Ey, Ez
         EN(emf),   E0(emf),   DE(emf),          // New E, old E, perturbation E
   
         dt(deltat),
         idx(static_cast< complex<double> >(0.5/deltax))
         {
// - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - 
            Nbc = Input::List().BoundaryCells;
            szx = Input::List().NxLocal[0];
    }

Here is the call graph for this function:

Here is the caller graph for this function:

Member Function Documentation

advance()

void Electric_Field_Methods::Implicit_E_Field::advance ( Algorithms::RK3< State1D > *  rk, State1D &  Y, collisions &  coll, VlasovFunctor1D_implicitE_p2 *  rkF  )

virtual

Implements Electric_Field_Methods::Efield_Method.

Definition at line 268 of file implicitE.cpp.

References collisions::advancef1(), Ampere(), Electric_Field_Methods::Current_xyz::calculate_J(), DE, Det33(), Detx33(), Dety33(), Detz33(), dt, State1D::EMF(), EN, Electric_Field_Methods::Efield_xyz::Ex(), EMF1D::Ex(), Electric_Field_Methods::Efield_xyz::Ey(), EMF1D::Ey(), Electric_Field_Methods::Efield_xyz::Ez(), EMF1D::Ez(), FindDE(), J0, J_Ex, J_Ey, J_Ez, JN, Electric_Field_Methods::Current_xyz::Jx(), Electric_Field_Methods::Current_xyz::Jy(), Electric_Field_Methods::Current_xyz::Jz(), and szx.

Referenced by Implicit_E_Field().

                                                                                                          {//, double time, double dt){
//--------------------------------------------------------------
//  Calculate the implicit electric field
//--------------------------------------------------------------

        int zeros_in_det(1);      // This counts the number of zeros in the determinant 
        int execution_attempt(0); // This counts the number of attempts to find invert the E-field
        State1D Yh(Yin);
        Yh = 0.0;

        // vector<SHarmonic1D> f00;
        // vector<SHarmonic1D> f10;vector<SHarmonic1D> f11;
        // vector<SHarmonic1D> f20;vector<SHarmonic1D> f21;vector<SHarmonic1D> f22;

        // for (size_t s(0); s < Yin.Species(); ++s){
        //     f00.push_back(Yin.SH(s,0,0));
        //     f10.push_back(Yin.SH(s,1,0));f11.push_back(Yin.SH(s,1,1));
        //     f20.push_back(Yin.SH(s,2,0));f21.push_back(Yin.SH(s,2,1));f22.push_back(Yin.SH(s,2,2));
        // }

        FindDE(Yin.EMF());                           //  Reset DE

// - - - - - - - - - - - - - - - - - - - - - -
        while ( (zeros_in_det > 0) && ( execution_attempt < 4) ) {  // Execute this loop at most twice
            zeros_in_det = 0;                                       // Count the zeros of the determinant
            ++execution_attempt;                                // Count the execusion attempts
// - - - - - - - - - - - - - - - - - - - - - -
// - - - - - - - - - - - - - - - - - - - - - -  
                                                        // Effect of E = 0 on f00, f10, f11
            
            coll.advancef1(Yin,Yh);                 // Collisions for f10, f11
            // for (size_t s(0); s < Yin.Species(); ++s){
            //     Yin.SH(s,1,0) = Yh.SH(s,1,0);
            //     Yin.SH(s,1,1) = Yh.SH(s,1,0);
            // }   
            J0.calculate_J(Yh);
            
            
            // for (size_t s(0); s < Yin.Species(); ++s){
            //     Yin.SH(s,0,0) = f00[s];
            //     Yin.SH(s,1,0) = f10[s];Yin.SH(s,1,1) = f11[s];
            // // Yin.SH(s,2,0) = f20[s];Yin.SH(s,2,1) = f21[s];Yin.SH(s,2,2) = f22[s];
            // }       

            Yin.EMF().Ex() = DE.Ex();             // Ex = DEx
            Yin = (*rk)(Yin, dt, rkF, 1);           // Effect of DEx on Y00, Y10, Y11, Y20, Y21, Y22            
            coll.advancef1(Yin,Yh);                 // Collisions for f10, f11
            // for (size_t s(0); s < Yin.Species(); ++s){
            //     Yin.SH(s,1,0) = Yh.SH(s,1,0);
            //     Yin.SH(s,1,1) = Yh.SH(s,1,0);
            // }                    // Collisions for f10, f11
            J_Ex.calculate_J(Yh);                // Evaluate J(DEx)
            
            // Ytemp = Yin;
            // for (size_t s(0); s < Yin.Species(); ++s){
            // Yin.SH(s,0,0) = f00[s];
            // Yin.SH(s,1,0) = f10[s];Yin.SH(s,1,1) = f11[s];
            // Yin.SH(s,2,0) = f20[s];Yin.SH(s,2,1) = f21[s];Yin.SH(s,2,2) = f22[s];
            // }       
           
                // 
            Yin.EMF().Ey() = DE.Ey();             // Ey = DEy
            Yin = (*rk)(Yin, dt, rkF, 2);           // Effect of DEy on Y00, Y10, Y11, Y20, Y21, Y22
            coll.advancef1(Yin,Yh);                 // Collisions for f10, f11
            // for (size_t s(0); s < Yin.Species(); ++s){
            //     Yin.SH(s,1,0) = Yh.SH(s,1,0);
            //     Yin.SH(s,1,1) = Yh.SH(s,1,0);
            // }                    // Collisions for f10, f11
            J_Ey.calculate_J(Yh);                 // Evaluate J(DEy)
            
            // Ytemp = Yin;
            // for (size_t s(0); s < Yin.Species(); ++s){
            //     Yin.SH(s,0,0) = f00[s];
            //     Yin.SH(s,1,0) = f10[s];Yin.SH(s,1,1) = f11[s];
            //     Yin.SH(s,2,0) = f20[s];Yin.SH(s,2,1) = f21[s];Yin.SH(s,2,2) = f22[s];
            // }       
// - - - - - - - - - - - - - - - - - - - - - -

           
            Yin.EMF().Ez() = DE.Ez();             // Ey = DEy
            Yin = (*rk)(Yin, dt, rkF, 3);           // Effect of DEy on Y00, Y10, Y11, Y20, Y21, Y22
            
            coll.advancef1(Yin,Yh);                 // Collisions for f10, f11
            // for (size_t s(0); s < Yin.Species(); ++s){
            //     Yin.SH(s,1,0) = Yh.SH(s,1,0);
            //     Yin.SH(s,1,1) = Yh.SH(s,1,0);
            // }   

            J_Ez.calculate_J(Yh); //exit(1);               // Evaluate J(DEy)
            
            // Ytemp = Yin;
            // for (size_t s(0); s < Yin.Species(); ++s){
            //     Yin.SH(s,0,0) = f00[s];
            //     Yin.SH(s,1,0) = f10[s];Yin.SH(s,1,1) = f11[s];
            //     Yin.SH(s,2,0) = f20[s];Yin.SH(s,2,1) = f21[s];Yin.SH(s,2,2) = f22[s];
            // }    
// - - - - - - - - - - - - - - - - - - - - - -

                
            Ampere(Yin.EMF());                           // Calculate JN
                
            
// - - - - - - - - - - - - - - - - - - - - - -
            
            //  calculate new EN 

            Array2D< complex<double> > sgm(3,3);     
            valarray< complex<double> > clm(3);

            
                // for (size_t iy(0); iy < szy; ++iy){
                    for (size_t ix(0); ix < szx; ++ix){
                        sgm(0,0) =  ( J_Ex.Jx()(ix) - J0.Jx()(ix) ) / DE.Ex()(ix);  // sxx = dJ(E_x)_x/DE_x
                        sgm(0,1) =  ( J_Ey.Jx()(ix) - J0.Jx()(ix) ) / DE.Ey()(ix);  // sxy = dJ(E_y)_x/DE_y
                        sgm(0,2) =  ( J_Ez.Jx()(ix) - J0.Jx()(ix) ) / DE.Ez()(ix);  // sxz = dJ(E_z)_x/DE_z

                        sgm(1,0) =  ( J_Ex.Jy()(ix) - J0.Jy()(ix) ) / DE.Ex()(ix);  // syx = dJ(E_x)_y/DE_x
                        sgm(1,1) =  ( J_Ey.Jy()(ix) - J0.Jy()(ix) ) / DE.Ey()(ix);  // syy = dJ(E_y)_y/DE_y
                        sgm(1,2) =  ( J_Ez.Jy()(ix) - J0.Jy()(ix) ) / DE.Ez()(ix);  // syz = dJ(E_z)_y/DE_z

                        sgm(2,0) =  ( J_Ex.Jz()(ix) - J0.Jz()(ix) ) / DE.Ex()(ix);  // szx = dJ(E_x)_z/DE_x
                        sgm(2,1) =  ( J_Ey.Jz()(ix) - J0.Jz()(ix) ) / DE.Ey()(ix);  // szy = dJ(E_y)_z/DE_y
                        sgm(2,2) =  ( J_Ez.Jz()(ix) - J0.Jz()(ix) ) / DE.Ez()(ix);  // szz = dJ(E_z)_z/DE_z
  
                        clm[0] =  JN.Jx()(ix) - J0.Jx()(ix);
                        clm[1] =  JN.Jy()(ix) - J0.Jy()(ix);
                        clm[2] =  JN.Jz()(ix) - J0.Jz()(ix);

                        

                        // std::cout << "\n\n" << sgm(0,0) << "  ,  " << sgm(0,1) << "  ,  " << sgm(0,2) << "\n";
                        // std::cout           << sgm(1,0) << "  ,  " << sgm(1,1) << "  ,  " << sgm(1,2) << "\n";
                        // std::cout           << sgm(2,0) << "  ,  " << sgm(2,1) << "  ,  " << sgm(2,2) << "\n";

                        
                    // Solve the 3 by 3 system of equations
                        complex<double> D_sgm( Det33(sgm) );            // The Determinant of the conductivity tensor
                        if ( abs( D_sgm.real() ) > 6.0*DBL_MIN ) {
                            EN.Ex()(ix) = Detx33(clm, sgm) / D_sgm;
                            EN.Ey()(ix) = Dety33(clm, sgm) / D_sgm;
                            EN.Ez()(ix) = Detz33(clm, sgm) / D_sgm;
                        }
                        else {
                            ++zeros_in_det;                                 // Try again with a new perturbation as below
                        // exit(1);
                            DE.Ex()(ix) *= 2;                        // DEx is 2 times larger than before
                            DE.Ey()(ix) *= 1.3;                      // DEy is 1.3 times larger  than before

                            EN.Ex()(ix) = 0.0; 
                            EN.Ey()(ix) = 0.0;
                            EN.Ez()(ix) = 0.0;
                        }
                    }
                // }
// - - - - - - - - - - - - - - - - - - - - - -

        } // End the while statement
// - - - - - - - - - - - - - - - - - - - - - -

// - - - - - - - - - - - - - - - - - - - - - -
        Yin.EMF().Ex() = EN.Ex();
        Yin.EMF().Ey() = EN.Ey();
        Yin.EMF().Ez() = EN.Ez();

        if ( zeros_in_det > 0) {
                        cout << "WARNING, Det = 0 in "<<zeros_in_det <<"locations" << endl;
                        if ( zeros_in_det > 8) { exit(1);}
        }

        // return Ytemp;
    }

Here is the call graph for this function:

Here is the caller graph for this function:

Ampere()

void Electric_Field_Methods::Implicit_E_Field::Ampere ( EMF1D &  emf )

private

Definition at line 484 of file implicitE.cpp.

References EMF1D::By(), EMF1D::Bz(), idx, JN, Electric_Field_Methods::Current_xyz::Jy(), and Electric_Field_Methods::Current_xyz::Jz().

Referenced by advance().

                                                                 {
//--------------------------------------------------------------
//  Use Ampere's law to calculate JN
//--------------------------------------------------------------
        Field1D tmpJi(emf.Bz());

//      Jx =  dBz/dy       
//      tmpJi    = Yin.EMF().Bz(); (same as assignment above)
        // tmpJi   *= idy;
        // JN.Jx()  = tmpJi.Dy();

//      Jy = -dBz/dx       
        tmpJi    = emf.Bz(); 
        tmpJi   *= (-1.0) * idx;
        tmpJi.Dx();
        // tmpJi(tmpJi.numx()-1) = 0.0;
        JN.Jy()  = tmpJi;        

//      Jz = -dBx/dy       
        // tmpJi    = Y.EMF().Bx(); 
        // tmpJi   *= (-1.0) * idy;
        // JN.Jz()  = tmpJi.Dy();    

//      Jz += dBy/dx       
        tmpJi    = emf.By(); 
        tmpJi   *=  idx;
        tmpJi.Dx();
        // tmpJi(tmpJi.numx()-1) = 0.0;
        JN.Jz() += tmpJi;   
    }

Here is the call graph for this function:

Here is the caller graph for this function:

FindDE()

void Electric_Field_Methods::Implicit_E_Field::FindDE ( EMF1D &  emf )

private

Definition at line 451 of file implicitE.cpp.

References DE, Electric_Field_Methods::Efield_xyz::Ex(), EMF1D::Ex(), Electric_Field_Methods::Efield_xyz::Ey(), EMF1D::Ey(), Electric_Field_Methods::Efield_xyz::Ez(), EMF1D::Ez(), and szx.

Referenced by advance().

                                                                 {
//--------------------------------------------------------------
//  Reset DE
//--------------------------------------------------------------
        double Eps(16.0*numeric_limits<double>::epsilon());  
        double LargeEps(sqrt(Eps));  

        DE.Ex() = emf.Ex();
        DE.Ey() = emf.Ey();
        DE.Ez() = emf.Ez();

        // Calculate DEx = LargeEps*(|E|)+Eps
        // for (size_t iy(0); iy < szy; ++iy){
            for (size_t ix(0); ix < szx; ++ix){
                DE.Ex()(ix) *=  DE.Ex()(ix);
                DE.Ey()(ix) *=  DE.Ey()(ix);
                DE.Ez()(ix) *=  DE.Ez()(ix);

                DE.Ex()(ix) +=  DE.Ey()(ix);
                DE.Ex()(ix) +=  DE.Ez()(ix);
                DE.Ex()(ix)  =  sqrt(DE.Ex()(ix)); 
                DE.Ex()(ix) *=  LargeEps;
                DE.Ex()(ix) +=  Eps;
            }
        // }
        DE.Ey() = DE.Ex(); 
        DE.Ez() = DE.Ex(); 

    }

Here is the call graph for this function:

Here is the caller graph for this function:

Field Documentation

DE

Efield_xyz Electric_Field_Methods::Implicit_E_Field::DE

private

Definition at line 165 of file implicitE.h.

Referenced by advance(), and FindDE().

dt

double Electric_Field_Methods::Implicit_E_Field::dt

private

Definition at line 166 of file implicitE.h.

Referenced by advance().

E0

Efield_xyz Electric_Field_Methods::Implicit_E_Field::E0

private

Definition at line 165 of file implicitE.h.

EN

Efield_xyz Electric_Field_Methods::Implicit_E_Field::EN

private

Definition at line 165 of file implicitE.h.

Referenced by advance().

idx

complex<double> Electric_Field_Methods::Implicit_E_Field::idx

private

Definition at line 167 of file implicitE.h.

Referenced by Ampere().

J0

Current_xyz Electric_Field_Methods::Implicit_E_Field::J0

private

Definition at line 164 of file implicitE.h.

Referenced by advance().

J_Ex

Current_xyz Electric_Field_Methods::Implicit_E_Field::J_Ex

private

Definition at line 164 of file implicitE.h.

Referenced by advance().

J_Ey

Current_xyz Electric_Field_Methods::Implicit_E_Field::J_Ey

private

Definition at line 164 of file implicitE.h.

Referenced by advance().

J_Ez

Current_xyz Electric_Field_Methods::Implicit_E_Field::J_Ez

private

Definition at line 164 of file implicitE.h.

Referenced by advance().

JN

Current_xyz Electric_Field_Methods::Implicit_E_Field::JN

private

Definition at line 164 of file implicitE.h.

Referenced by advance(), and Ampere().

Nbc

int Electric_Field_Methods::Implicit_E_Field::Nbc

private

Definition at line 184 of file implicitE.h.

Referenced by Implicit_E_Field().

szx

int Electric_Field_Methods::Implicit_E_Field::szx

private

Definition at line 184 of file implicitE.h.

Referenced by advance(), FindDE(), and Implicit_E_Field().

---

The documentation for this class was generated from the following files:

• /Users/archis/Dropbox/work/dev/oshun/oshun-OS/OSHUN/1d_cpp/source/implicitE.h
• /Users/archis/Dropbox/work/dev/oshun/oshun-OS/OSHUN/1d_cpp/source/implicitE.cpp