NewProcess

Revision as of 18:30, 28 Nov 2005
217.20.209.178 (Talk | contribs)

← Go to previous diff
Revision as of 20:40, 21 Jan 2006
Hatta (Talk | contribs)

Go to next diff →
Line 14: Line 14:
Для задействования внешних сечений было решено написать свой процесс, поскольку тольщина таллиевого образца была достаточно солидной, чтобы гамма-кванты порядочно застревали в нем до того, как произойдет фотоядерная реакция. Для задействования внешних сечений было решено написать свой процесс, поскольку тольщина таллиевого образца была достаточно солидной, чтобы гамма-кванты порядочно застревали в нем до того, как произойдет фотоядерная реакция.
 +
 +==G4DiscreteProcess==
 +
 +В моем случае можно было пренебречь вторичными фотонами, которые могут родиться в результате фотоядерной реакции, т. к. вероятность для такого фотона вновь провзаимодействовать с ядром ничтожна, а сами по себе фотоны меня не интересовали. Кроме того, для корректного их учета необходимо как минимум знать угловые сечения для данных реакций, а у меня было только интегральное сечение вылета одного или нескольких нейтронов — узнать направление вторичного фотона по этим данным мне представляется затруднительным. То же самое можно сказать и о нейтронах: я не знаю куда они полетят и с какой энергией, и, кроме того, мне это не важно. Поэтому решено было просто учитывать факт реакции, а затем просто убивать исходный фотон.
 +
 +Ядерная реакция подходит под критерий дискретного процесса (см. в документации GEANT про типы процессов), поэтому для реализации необходимо создать класс, производный от <tt>G4VDiscreteProcess</tt>.
 +
 + class CuteProcess: public G4VDiscreteProcess
 + {
 + public:
 + /// Конструктор.
 + /// @param dirname: директория с сечениями;
 + /// @param nn: число выделяемых нейтронов;
 + CuteProcess(const G4String& name,
 + const G4String& dirname,
 + int nn):
 + G4VDiscreteProcess(name, fUserDefined),
 + nneutrons(nn)
 + {
 + loadFromDir(dirname, nn);
 + }
 + virtual G4VParticleChange *PostStepDoIt(const G4Track& track, const G4Step& step);
 +protected:
 + virtual G4double GetMeanFreePath(const G4Track& aTrack,
 + G4double previousStepSize,
 + G4ForceCondition* condition);
 + void loadFromDir(const G4String& dirname, int nn);
 +private:
 + Maths::Interpolation::Linear* xs;
 + int nneutrons;
 + };
 +
 +Как видно, в конструкторе не делается ничего особенно важного, только загружаются сами сечения из файлов. Я не буду вдаваться в детали интерполяции значений сечений, их формата на диске и т.&nbsp;д., а остановлюсь на действительно важных деталях. Два метода в этом классе, которые и реализуют всю физику процесса это <tt>PostStepDoIt</tt> и <tt>GetMeanFreePath</tt>. Последний вызывается каждый раз при очередном шажке частицы в симуляции и определяет вероятность вызова данного конкретного процесса. <tt>PostStepDoIt</tt>&nbsp;&mdash; это то, что делает процесс, когда он все-таки происходит.
 +
 +
 + G4double CuteProcess::GetMeanFreePath(const G4Track& aTrack,
 + G4double previousStepSize,
 + G4ForceCondition* condition)
 + {
 + double energy = aTrack.GetKineticEnergy();
 + if (aTrack.GetMaterial()->GetName() != "Thallium"
 + || aTrack.GetDefinition() != G4Gamma::GammaDefinition())
 + return 1e100;
 + double crosssec = xs->getValue(energy); // Сечение из файла.
 + double l = 1/(6.022e23/mole*11.85*g/cm3/(204.3833*g/mole)*crosssec);
 + if(l<0) l = 1e100;
 + return l;
 + }
 +
 +Здесь по простой формуле вычисляется из сечения свободный пробег, который и возвращается. <tt>crossec</tt> уже переведен во внутренние единицы GEANT4 (кв. мм). Проверка с if на второй строчке это всего лишь костыль, спасающий от реальности: даны сечения только для таллия, так что ни в каких других материалах процесс не будет активизирован, т.&nbsp;к. в этом случае он вернет пробег 1e100.
 +
 +Вот, что мы делаем, когда процесс-таки вызывается:
 +
 + G4VParticleChange* CuteProcess::PostStepDoIt(const G4Track& track, const G4Step& step)
 + {
 + aParticleChange.Initialize(track);
 + aParticleChange.ProposeEnergy(0.);
 + aParticleChange.ProposeTrackStatus(fStopAndKill);
 + std::cout << "gen_neutron: " << nneutrons << " more" << std::endl;
 + return &aParticleChange;
 + }
 +
 +<tt>aParticleChange</tt>&nbsp;&mdash; это член класса '''G4VDiscreteProcess''', в котором вы должны описать все изменения, которые вызываются вашим процессом. В данном случае, как видно, я устанавливаю статус трекинга частицы в <tt>fStopAndKill</tt>, что соответствует выводу ее из процесса симуляции. Установка энергии строчкой выше, на самом деле не нужна, и здесь только ради перестраховки. Кроме того, я вывожу в '''stdout''' сообщение о срабатывании процесса.
 +
 +Подключается такой самодельный процесс точно так же, как и встроенные:
 +
 + pmanager->AddDiscreteProcess(new CuteProcess("AyAy", "/home/hatta/work/203Tl/", n));

Revision as of 20:40, 21 Jan 2006

Введение

Предупреждаю! Эта статья еще не дописана.


Geant знает достаточно большое количество физических моделей. Они в него зашиты, как, просто, основная часть кода. Зайдите интереса ради в папку $G4INSTALL/source/processes и оглядитесь вокруг&nbsp— одни только названия файлов о многом расскажут. Но, несмотря на обилие встроенных моделей, иногда может понадобиться примитивная в сущности вещь: сделать так, чтобы Geant считал по вашим сечениям. Об том, как это сделать, и написана эта статья.

Проверено на себе

Лично мне пришлось заняться этим на первый взгляд бессмысленным делом, когда руководство разбраковало в пух и прах встроенные в Geant адронные модели и дало мне несколько файлов с фразой «Вот тебе надежные сечения». Не объяснив при этом даже формат хранения данных… но это к делу не относится.

Поглядев в файлах, я установил, что сечения были даны в виде таблицы «Энергия—Сечение». Всего было 6 файлов для 6 реакций: (gamma, n), (gamma, 2n) и так далее. Да, и все это было для одного единственного ядра, 203Tl (таллий).

Для задействования внешних сечений было решено написать свой процесс, поскольку тольщина таллиевого образца была достаточно солидной, чтобы гамма-кванты порядочно застревали в нем до того, как произойдет фотоядерная реакция.

G4DiscreteProcess

В моем случае можно было пренебречь вторичными фотонами, которые могут родиться в результате фотоядерной реакции, т. к. вероятность для такого фотона вновь провзаимодействовать с ядром ничтожна, а сами по себе фотоны меня не интересовали. Кроме того, для корректного их учета необходимо как минимум знать угловые сечения для данных реакций, а у меня было только интегральное сечение вылета одного или нескольких нейтронов — узнать направление вторичного фотона по этим данным мне представляется затруднительным. То же самое можно сказать и о нейтронах: я не знаю куда они полетят и с какой энергией, и, кроме того, мне это не важно. Поэтому решено было просто учитывать факт реакции, а затем просто убивать исходный фотон.

Ядерная реакция подходит под критерий дискретного процесса (см. в документации GEANT про типы процессов), поэтому для реализации необходимо создать класс, производный от G4VDiscreteProcess.

class CuteProcess: public G4VDiscreteProcess
{
public:
 /// Конструктор.
 /// @param dirname: директория с сечениями;
 /// @param nn: число выделяемых нейтронов;
 CuteProcess(const G4String& name,
          const G4String& dirname,
          int nn): 
      G4VDiscreteProcess(name, fUserDefined),
      nneutrons(nn)
 {
   loadFromDir(dirname, nn);
 }
 virtual G4VParticleChange *PostStepDoIt(const G4Track& track, const G4Step& step);

protected:

 virtual G4double GetMeanFreePath(const G4Track& aTrack,
                                  G4double previousStepSize,
                                  G4ForceCondition* condition);
 void loadFromDir(const G4String& dirname, int nn);

private:

 Maths::Interpolation::Linear* xs;
 int nneutrons;
};

Как видно, в конструкторе не делается ничего особенно важного, только загружаются сами сечения из файлов. Я не буду вдаваться в детали интерполяции значений сечений, их формата на диске и т. д., а остановлюсь на действительно важных деталях. Два метода в этом классе, которые и реализуют всю физику процесса это PostStepDoIt и GetMeanFreePath. Последний вызывается каждый раз при очередном шажке частицы в симуляции и определяет вероятность вызова данного конкретного процесса. PostStepDoIt — это то, что делает процесс, когда он все-таки происходит.


G4double CuteProcess::GetMeanFreePath(const G4Track& aTrack,
                                  G4double previousStepSize,
                                  G4ForceCondition* condition)
{
 double energy = aTrack.GetKineticEnergy();
 if (aTrack.GetMaterial()->GetName() != "Thallium" 
  || aTrack.GetDefinition() != G4Gamma::GammaDefinition())
   return 1e100;
 double crosssec = xs->getValue(energy); // Сечение из файла.
 double l = 1/(6.022e23/mole*11.85*g/cm3/(204.3833*g/mole)*crosssec);
 if(l<0) l = 1e100;
 return l;
}

Здесь по простой формуле вычисляется из сечения свободный пробег, который и возвращается. crossec уже переведен во внутренние единицы GEANT4 (кв. мм). Проверка с if на второй строчке это всего лишь костыль, спасающий от реальности: даны сечения только для таллия, так что ни в каких других материалах процесс не будет активизирован, т. к. в этом случае он вернет пробег 1e100.

Вот, что мы делаем, когда процесс-таки вызывается:

G4VParticleChange* CuteProcess::PostStepDoIt(const G4Track& track, const G4Step& step)
{
 aParticleChange.Initialize(track);
 aParticleChange.ProposeEnergy(0.);
 aParticleChange.ProposeTrackStatus(fStopAndKill);
 std::cout << "gen_neutron: " << nneutrons << " more" << std::endl;
 return &aParticleChange;
}

aParticleChange — это член класса G4VDiscreteProcess, в котором вы должны описать все изменения, которые вызываются вашим процессом. В данном случае, как видно, я устанавливаю статус трекинга частицы в fStopAndKill, что соответствует выводу ее из процесса симуляции. Установка энергии строчкой выше, на самом деле не нужна, и здесь только ради перестраховки. Кроме того, я вывожу в stdout сообщение о срабатывании процесса.

Подключается такой самодельный процесс точно так же, как и встроенные:

pmanager->AddDiscreteProcess(new CuteProcess("AyAy", "/home/hatta/work/203Tl/", n));
Edit page