spatialg.cpp

pic 模拟程序！面向对象

  
  region->setParticleLimit( CP.get_particleLimit() );
  region->set_np2cFactor( CP.get_np2cFactor() );
  region->setduplicateParticles( (int)(CP.getduplicateParticles()+0.1) );

  region->set_histmax(CP.get_histmax());
  
#ifdef MPI_VERSION
  char buf[512];
  if(Name.getValue()==ostring("Noname") )
    sprintf(buf,"Region %d",MPI_RANK);
  else
    sprintf(buf,"%s Region %d",Name.getValue().c_str(),MPI_RANK);
  region->setName(buf);
#else  /* MPI_VERSION */
  region->setName(Name.getValue());
#endif  /* MPI_VERSION */
  oopicList<Diag> *DList = new DiagList;
  oopicListIter<DiagParams> nD(Diags);
  for (nD.restart(); !nD.Done(); nD++){
    try{
      nD.current()->toGrid(G);
    }
    catch(Oops& oops){
      oops.prepend("SpatialRegionGroup::CreateCounterPart: Error:\n"); //OK
      throw oops;
    }

#ifdef MPI_VERSION
    if(! (nD.current()->j1.getValue()==-1 && nD.current()->A1.getValue()==-1)) {
      // we don't want to create boundaries which are out-of-bounds!
#endif /* MPI_VERSION */
      Diag* dPtr = 0;
      try{
        dPtr = nD.current()->CreateCounterPart(region);
      }
      catch(Oops& oops2){
        oops2.prepend("SpatialRegionGroup::CreateCounterPart: Error: \n");//OK
        throw oops2;
      }
      DList->add(dPtr);
#ifdef MPI_VERSION
    }
#endif
    
  }
  region->setDiagList(DList);

// create a pointer to list of H5Diagnostic parameters

  oopicList<H5DiagParams> *H5DList = new oopicList<H5DiagParams>;
  oopicListIter<H5DiagParams> nDD(H5DiagDumpParams);
  for (nDD.restart(); !nDD.Done(); nDD++){
     H5DList->add(nDD.current());
  }
  region->setH5DiagDumpList(H5DList);
    
  if (WasDumpFileGiven)	{
    try{
      Restore(region);
    }
    catch(Oops& oops2){
      oops2.prepend("SpatialRegionGroup::CreateCounterPart: Error: \n");//OK
      throw oops2;
    }
  }
  
  //	else 
  try {
    ApplyToList(CreateCounterPart(region), Loads, LoadParams);
  }
  catch(Oops& oops2){
      oops2.prepend("SpatialRegion::CreateCounterPart: Error: \n");
      throw oops2;
  }
  
  F->compute_rho(region->getParticleGroupListArray(),
                 region->get_nSpecies());
  
#ifdef MPI_VERSION
  linkSRB();
#endif
  try{
    ApplyToList(putCharge_and_Current(region->getTime()),*BList,Boundary);
  }
  catch(Oops& oops2){
      oops2.prepend("SpatialRegion::CreateCounterPart: Error: \n");//done
      throw oops2;
  }
  
  // Need to receive a message sent by putCharge_and_Current
  // in MPI, to set up toNodes applyFields with dt=0 just receives
  // sent fields and sets up the SRB for toNodes
#ifdef MPI_VERSION
  try{
    ApplyToList(applyFields(region->getTime(),0),SpatialRegionBoundaryList, SpatialRegionBoundary);
  }
  catch(Oops& oops3){
    oops3.prepend("SpatialRegionGroup::CreateCounterPart: Error: \n");
    throw oops3;
  }

#endif

   // Some boundaries need feedback info from SpatialRegion
   //
   //
   ApplyToList(setRegion(region),*BList,Boundary);
   
   /* this line need to follow compute_rho so that the boundary 
      get rho correct */

   if ((CP.getElectrostaticFlag() || !WasDumpFileGiven) && CP.getInitPoissonSolveFlag()){

#ifndef MPI_VERSION  /* ES field not working in MPI version yet! */
    F->ElectrostaticSolve(region->getTime(),region->get_dt());
#endif /* MPI_VERSION */
    
  }

#ifdef MPI_VERSION /* Since the ES field solve is not working, at least apply BC */
  if( CP.getInitPoissonSolveFlag()) {
    try{
      ApplyToList(applyFields(region->getTime(),region->get_dt()),*BList,Boundary);
    }
    catch(Oops& oops3){
      oops3.prepend("SpatialRegionGroup::CreateCounterPart: Error: \n");
      throw oops3;
    }
  }
#endif
  
  F->setdt(region->get_dt());
  
  //  the next two lines set up initial EM fields
  //  which are added to any previous fields.
  //  since it is additive, we should do it only once.
  if(!WasDumpFileGiven) {
    F->setEinit(CP.getE1init(),CP.getE2init(),CP.getE3init());
    F->setBinit(CP.getB1init(),CP.getB2init(),CP.getB3init());
  }
  
  F->toNodes(region->getTime());
  
  
  // set up the moving window.
  switch(CP.getMovingWindow()){
  case 1:
    region->setShift(Vector2(-1,0));
    break;
  case 2:
    region->setShift(Vector2(1,0));
    break;
  case 3:
    region->setShift(Vector2(-1,0));
    break;
  case 4:
    region->setShift(Vector2(1,0));
    break;
  }

  region->setShiftDelayTime(CP.getShiftDelayTime());
  region->setGasOffTime(CP.getGasOffTime());

  
  return region;
}


//-------------------------------------------------------------------
//	temp mechanism to show rule messages - should work on the entire
//	SRG?  Should SRG inherit from ParameterGroup?

void	SpatialRegionGroup::showRuleMessages()
{
   if (GP.ruleMessages.nitems) GP.showRuleMessages();
   if (CP.ruleMessages.nitems) CP.showRuleMessages();
   oopicListIter<BoundaryParams> nB(Boundaries);
   for (nB.restart(); !nB.Done(); nB++)
      {
         BoundaryParams *bp = nB.current();
         if (bp->ruleMessages.nitems) bp->showRuleMessages();
      }
}                               

// the following are used in GeometryEditor::setSpatialRegion
//
GridParams* SpatialRegionGroup::getGridParams()
{
   return &GP;
}

ControlParams* SpatialRegionGroup::getControlParams()
{
   return &CP;
}

oopicList<BoundaryParams>* SpatialRegionGroup::getBoundaryParamsList()
{
   return &Boundaries;
}

oopicList<DiagnosticParams>* SpatialRegionGroup::getProbeParamsList()
{
   return &Probes;
}

// the following are used in GeometryEditor createDragObject
BoundaryParams* SpatialRegionGroup::constructBoundaryParams(ostring groupname)
  throw(Oops){
   BoundaryParams* bptr = NULL;

   if (groupname == (ostring)"Conductor")
      bptr = new ConductorParams(&GP,this);

   if (groupname == (ostring)"Polarizer")
      bptr = new PolarParams(&GP,this);

   if (groupname == (ostring)"SpatialRegionBoundary"){
     try{    
       bptr = new SpatialRegionBoundaryParams(&GP,this);
     }
     catch(Oops& oops){
       oops.prepend("SpatialRegionGroup::constructBoundaryParams: Error:\n"); //SpatialRegionGroup::CreateBoundaryParams
       throw oops;
     }
   }

   if (groupname == (ostring)"CylindricalAxis")
      bptr = new CylindricalAxisParams(&GP,this);

   if (groupname == (ostring)"FieldEmitter")
      bptr = new FieldEmitterParams(&GP,this);

   if (groupname == (ostring)"FieldEmitter2")
      bptr = new FieldEmitterParams2(&GP,this);

   if (groupname == (ostring)"BeamEmitter")
      bptr = new BeamEmitterParams(&GP,this);

   if (groupname == (ostring)"FowlerNordheimEmitter")
      bptr = new FowlerNordheimEmitterParams(&GP,this);

   if (groupname == (ostring)"SingleParticleEmitter")
      bptr = new SingleParticleEmitterParams(&GP,this);

   if (groupname == (ostring)"PlasmaSource")
      bptr = new PlasmaSourceParams(&GP,this);

   if (groupname == (ostring)"LaserPlasmaSource")
      bptr = new LaserPlasmaSourceParams(&GP,this);

   if (groupname == (ostring)"EmitPort")
      bptr = new EmitPortParams(&GP,this);
   
   if (groupname == (ostring)"VarWeightBeamEmitter")
      bptr = new VarWeightBeamEmitterParams(&GP,this);

   if (groupname == (ostring)"ExitPort")
      bptr = new ExitPortParams(&CP, &GP,this);

   if (groupname == (ostring)"Iloop")
      bptr = new IloopParams(&GP, this);

   if (groupname == (ostring)"CurrentRegion")
      bptr = new CurrentRegionParams(&GP,this);

   if (groupname == (ostring)"PortTE")
      bptr = new PortTEParams(&GP,this);

   if (groupname == (ostring)"PortTM")
      bptr = new PortTMParams(&GP,this);

   if (groupname == (ostring)"PortTEM_MW")
      bptr = new PortTEM_MW_Params(&GP,this);

   if (groupname == (ostring)"PortGauss")
      bptr = new PortGaussParams(&GP,this);

   if (groupname == (ostring)"AbsorbWave")
      bptr = new AbsorbWaveParams(&GP,this);

   if (groupname == (ostring)"LaunchWave")
      bptr = new LaunchWaveParams(&GP,this);

   if (groupname == (ostring)"Gap")
      bptr =	new GapParams(&GP,this);

   if (groupname == (ostring)"Equipotential")
      bptr = new EquipotentialParams(&GP,this);
   
   if (groupname == (ostring)"Dielectric")
      bptr = new DielectricParams(&GP,this);

   if (groupname == (ostring)"DielectricRegion")
      bptr = new DielectricRegionParams(&GP,this);

   if (groupname == (ostring)"DielectricTriangle")
      bptr = new DielectricTriangleParams(&GP,this);

   if (groupname == (ostring)"Foil")
      bptr = new FoilParams(&GP,this);

   if (groupname == (ostring)"CurrentSource")
     bptr = new CurrentSourceParams(&GP,this);

   if (bptr) {
     bptr->setName(groupname);
   }
   else {
     cerr << "\n\nWARNING:  Invalid boundary type!" << endl;
     cerr << "  " << groupname << " does not exist..." << endl;
     cerr << "  You must modify your input file." << endl;
   }
   return bptr;
}

DiagnosticParams* SpatialRegionGroup::constructProbeParams(ostring groupname)
{
   DiagnosticParams* dptr=0;
#ifndef TEXT_OOPIC
   /*
      #ifndef UNIX
      if (groupname == (ostring)"E Field")
      dptr = new EFieldProbe(&GP);
      if (groupname == (ostring)"B Field")
      dptr = new BFieldProbe(&GP);
      if (groupname == (ostring)"Current")
      dptr = new CurrentProbe(&GP);
      if (groupname == (ostring)"DivDerror")
      dptr = new DivDerrorProbe(&GP);
      dptr->setName(groupname);
      #endif
      */
#endif
   return dptr;
}

ostring 	SpatialRegionGroup::readSpecies(ifstream& fin) 
{
   return advisorp->readSpecies(fin);
}

BOOL SpatialRegionGroup::KnownSpecies(ostring &speciesName)
{
   return advisorp->KnownSpecies(speciesName);
}

Species* SpatialRegionGroup::get_speciesPtr(const ostring &speciesName)
{
   return advisorp->get_speciesPtr(speciesName);
}

SpeciesList* SpatialRegionGroup::get_speciesList()
{
   return advisorp->get_speciesList();
}

SpatialRegionGroup::~SpatialRegionGroup() {
   Boundaries.deleteAll();
   //	Loads.deleteAll();   for moving window, want to keep these around.
   Diags.deleteAll();
   Probes.deleteAll();
   MCCParamsList.deleteAll();
};

#ifdef MPI_VERSION

void SpatialRegionGroup::addNeededSRBstoRegion(oopicList<Boundary> *BList,Grid *G) 
throw (Oops){
   // Grid has a mask in it so
   // I can figure out where the SRB's need to go!
   int J = G->getJ();
   int K = G->getK();
   int j,k;
   int bcount;  // count of SRB along an edge, needed for naming
   int Js,Ks;   // beginning of the current SRB
   int makingSRB=0;  // flag on whether we're making an SRB or not.
   //  0->1 transitions set Js, Ks,
   //  1->0 transitions create the boundary.
   int Nindex;  // index to neighbor array
   int normal;  // normal
   Boundary ***bc1=G->GetBC1();
   Boundary ***bc2=G->GetBC2();
   Boundary *newSRB;