Message: Help setting up an Electric Field Not Logged In (login)
 Next-in-Thread Next-in-Thread
 Next-in-Forum Next-in-Forum

None Help setting up an Electric Field 

Forum: Fields: Magnetic and Otherwise
Date: 10 Jul, 2007
From: Adam <ax_blais@laurentian.ca>

I have a xenon cylinder surrounded by air. I want to make an electric field in only the xenon. I tried borrowing the class from the field02 extended example. My field setup class is this:

#include "ExN02ElectricFieldSetup.hh"
#include "ExN02FieldMessenger.hh"

#include "G4UniformElectricField.hh"
#include "G4UniformMagField.hh"
#include "G4MagneticField.hh"
#include "G4FieldManager.hh"
#include "G4TransportationManager.hh"
#include "G4EquationOfMotion.hh"
#include "G4EqMagElectricField.hh"
#include "G4Mag_UsualEqRhs.hh"
#include "G4MagIntegratorStepper.hh"
#include "G4MagIntegratorDriver.hh"
#include "G4ChordFinder.hh"

#include "G4ExplicitEuler.hh"
#include "G4ImplicitEuler.hh"
#include "G4SimpleRunge.hh"
#include "G4SimpleHeum.hh"
#include "G4ClassicalRK4.hh"
#include "G4HelixExplicitEuler.hh"
#include "G4HelixImplicitEuler.hh"
#include "G4HelixSimpleRunge.hh"
#include "G4CashKarpRKF45.hh"
#include "G4RKG3_Stepper.hh"

////////////////////////////////////////////////////////////////////////// // // Constructors:

ExN02ElectricFieldSetup::ExN02ElectricFieldSetup()
  : fChordFinder(0), fStepper(0), fIntgrDriver(0)
{ 
  fEMfield = new G4UniformElectricField(
                 G4ThreeVector(0.0,1000.0*kilovolt/cm,0.0));
  fFieldMessenger = new ExN02FieldMessenger(this) ;
  fEquation = new G4EqMagElectricField(fEMfield); 
  fMinStep     = 0.010*mm ; // minimal step of 10 microns
  fStepperType = 4 ;        // ClassicalRK4 -- the default stepper

  fFieldManager = GetGlobalFieldManager();
  UpdateField();

}

/////////////////////////////////////////////////////////////////////////////////
ExN02ElectricFieldSetup::ExN02ElectricFieldSetup(G4ThreeVector pFV)
  : fChordFinder(0),
    fStepper(0),
    fIntgrDriver(0)
{
  fEMfield = new G4UniformElectricField(pFV);
  // GetGlobalFieldManager()->CreateChordFinder(this);

  fFieldMessenger = new ExN02FieldMessenger(this) ;
  fEquation = new G4EqMagElectricField(fEMfield); 
  fMinStep     = 0.010*mm ; // minimal step of 10 microns
  fStepperType = 4 ;        // ClassicalRK4 -- the default stepper

  fFieldManager = GetGlobalFieldManager();
  UpdateField();
}

////////////////////////////////////////////////////////////////////////////////

ExN02ElectricFieldSetup::~ExN02ElectricFieldSetup()
{ 
  if(fChordFinder) delete fChordFinder;
  if(fStepper)     delete fStepper;
  if(fEquation)    delete fEquation;
  if(fEMfield)     delete fEMfield;
}

///////////////////////////////////////////////////////////////////////////// // // Register this field to 'global' Field Manager and // Create Stepper and Chord Finder with predefined type, minstep (resp.) //

void ExN02ElectricFieldSetup::UpdateField()
{ 
  SetStepper();

  G4cout<<"The minimal step is equal to "<<fMinStep/mm<<" mm"<<G4endl ;

  fFieldManager->SetDetectorField(fEMfield );

  if(fChordFinder) delete fChordFinder;
  // fChordFinder = new G4ChordFinder( fEMfield, fMinStep, fStepper);

  fIntgrDriver = new G4MagInt_Driver(fMinStep,
                                     fStepper, 
                                     fStepper->GetNumberOfVariables() );

  fChordFinder = new G4ChordFinder(fIntgrDriver);

fFieldManager->SetChordFinder( fChordFinder ); }

///////////////////////////////////////////////////////////////////////////// // // Set stepper according to the stepper type //

void ExN02ElectricFieldSetup::SetStepper()
{ 
  G4int nvar = 8;

  if(fStepper) delete fStepper;

  switch ( fStepperType )
  { 
    case 0:  
      fStepper = new G4ExplicitEuler( fEquation, nvar );
      G4cout<<"G4ExplicitEuler is calledS"<<G4endl;
      break;
    case 1:  
      fStepper = new G4ImplicitEuler( fEquation, nvar );
      G4cout<<"G4ImplicitEuler is called"<<G4endl;
      break;
    case 2:  
      fStepper = new G4SimpleRunge( fEquation, nvar );
      G4cout<<"G4SimpleRunge is called"<<G4endl;
      break;
    case 3:  
      fStepper = new G4SimpleHeum( fEquation, nvar );
      G4cout<<"G4SimpleHeum is called"<<G4endl;
      break;
    case 4:  
      fStepper = new G4ClassicalRK4( fEquation, nvar );    
      G4cout<<"G4ClassicalRK4 (default) is called"<<G4endl;
      break;
    case 5:  
      fStepper = new G4CashKarpRKF45( fEquation, nvar );
      G4cout<<"G4CashKarpRKF45 is called"<<G4endl;
      break;
    case 6:  
      fStepper = 0; // new G4RKG3_Stepper( fEquation, nvar );       
      G4cout<<"G4RKG3_Stepper is not currently working for Electric Field"<<G4endl;     
      break;
    case 7:  
      fStepper = 0; // new G4HelixExplicitEuler( fEquation ); 
      G4cout<<"G4HelixExplicitEuler is not valid for Electric Field"<<G4endl;     
      break;
    case 8:  
      fStepper = 0; // new G4HelixImplicitEuler( fEquation ); 
      G4cout<<"G4HelixImplicitEuler is not valid for Electric Field"<<G4endl;     
      break;
    case 9:  
      fStepper = 0; // new G4HelixSimpleRunge( fEquation );   
      G4cout<<"G4HelixSimpleRunge is not valid for Electric Field"<<G4endl;
      break; 
    default: fStepper = 0;
  }
}

///////////////////////////////////////////////////////////////////////////// // // Set the value of the Global Field to fieldValue along Z //

void ExN02ElectricFieldSetup::SetFieldValue(G4double fieldValue)
{ 
  G4ThreeVector fieldVector( 0.0, 0.0, fieldValue );

SetFieldValue( fieldVector ); }

/////////////////////////////////////////////////////////////////////////////// // // Set the value of the Global Field value to fieldVector //

void ExN02ElectricFieldSetup::SetFieldValue(G4ThreeVector fieldVector)
{ 
  // Find the Field Manager for the global field
  G4FieldManager* fieldMgr= GetGlobalFieldManager();

  if(fieldVector != G4ThreeVector(0.,0.,0.))
  { 
    if(fEMfield) delete fEMfield;
    fEMfield = new  G4UniformElectricField(fieldVector);

    fEquation->SetFieldObj(fEMfield);  // must now point to the new field

    // UpdateField();

    fieldMgr->SetDetectorField(fEMfield);
  }
  else
  { 
    // If the new field's value is Zero, then it is best to
    //  insure that it is not used for propagation.
    if(fEMfield) delete fEMfield;
    fEMfield = 0;
    fEquation->SetFieldObj(fEMfield);   // As a double check ...

    G4MagneticField* fEMfield = 0;
    fieldMgr->SetDetectorField(fEMfield);
  }
}
////////////////////////////////////////////////////////////////////////////////
//
//  Utility method

G4FieldManager*  ExN02ElectricFieldSetup::GetGlobalFieldManager()
{
  return G4TransportationManager::GetTransportationManager()
         ->GetFieldManager();
}

I added the line:

ExN02ElectricFieldSetup* field = new ExN02ElectricFieldSetup();

to my main. When the simulation starts, it seems as though particles are leaving the xenon and keep going on forever until they exit the world volume, which gives me a segmentation fault. I know that to isolate the field to only the xenon cylinder, I need to use the line:

xenon_log->SetFieldManager(field->GetGlobalFieldManager(), false);

Is this right? If it is right, in what class do I put this line? The documentation isn't very specific about this.

Inline Depth:
 1 1
 All All
Outline Depth:
 1 1
 2 2
 All All
Add message: (add)

1 Feedback: Re: Help setting up an Electric Field   (Peter Gumplinger - 10 Jul, 2007)
(_ None: Re: Help setting up an Electric Field   (Adam - 10 Jul, 2007)
(_ None: Re: Help setting up an Electric Field   (John Apostolakis - 10 Jul, 2007)
(_ None: Re: Help setting up an Electric Field   (Adam - 11 Jul, 2007)
1 Warning: Re: Help setting up an Electric Field   (Peter Gumplinger - 11 Jul, 2007)
3 None: Re: Help setting up an Electric Field   (John Apostolakis - 13 Jul, 2007)
1 None: Re: Help setting up an Electric Field   (Adam - 13 Jul, 2007)
...
 Add Message Add Message
to: "Help setting up an Electric Field"

 Subscribe Subscribe

This site runs SLAC HyperNews version 1.11-slac-98, derived from the original HyperNews


[ Geant 4 Home | Geant 4 HyperNews | Search | Request New Forum | Feedback ]