| Geant4 Cross Reference |
1 // 1 //
2 // ******************************************* 2 // ********************************************************************
3 // * License and Disclaimer 3 // * License and Disclaimer *
4 // * 4 // * *
5 // * The Geant4 software is copyright of th 5 // * The Geant4 software is copyright of the Copyright Holders of *
6 // * the Geant4 Collaboration. It is provided 6 // * the Geant4 Collaboration. It is provided under the terms and *
7 // * conditions of the Geant4 Software License 7 // * conditions of the Geant4 Software License, included in the file *
8 // * LICENSE and available at http://cern.ch/ 8 // * LICENSE and available at http://cern.ch/geant4/license . These *
9 // * include a list of copyright holders. 9 // * include a list of copyright holders. *
10 // * 10 // * *
11 // * Neither the authors of this software syst 11 // * Neither the authors of this software system, nor their employing *
12 // * institutes,nor the agencies providing fin 12 // * institutes,nor the agencies providing financial support for this *
13 // * work make any representation or warran 13 // * work make any representation or warranty, express or implied, *
14 // * regarding this software system or assum 14 // * regarding this software system or assume any liability for its *
15 // * use. Please see the license in the file 15 // * use. Please see the license in the file LICENSE and URL above *
16 // * for the full disclaimer and the limitatio 16 // * for the full disclaimer and the limitation of liability. *
17 // * 17 // * *
18 // * This code implementation is the result 18 // * This code implementation is the result of the scientific and *
19 // * technical work of the GEANT4 collaboratio 19 // * technical work of the GEANT4 collaboration. *
20 // * By using, copying, modifying or distri 20 // * By using, copying, modifying or distributing the software (or *
21 // * any work based on the software) you ag 21 // * any work based on the software) you agree to acknowledge its *
22 // * use in resulting scientific publicati 22 // * use in resulting scientific publications, and indicate your *
23 // * acceptance of all terms of the Geant4 Sof 23 // * acceptance of all terms of the Geant4 Software license. *
24 // ******************************************* 24 // ********************************************************************
25 // 25 //
26 // 26 //
27 // Author: Mathieu Karamitros (kara (AT) cenbg 27 // Author: Mathieu Karamitros (kara (AT) cenbg . in2p3 . fr)
28 // 28 //
29 // WARNING : This class is released as a proto 29 // WARNING : This class is released as a prototype.
30 // It might strongly evolve or even disapear i 30 // It might strongly evolve or even disapear in the next releases.
31 // 31 //
32 // History: 32 // History:
33 // ----------- 33 // -----------
34 // 10 Oct 2011 M.Karamitros created 34 // 10 Oct 2011 M.Karamitros created
35 // 35 //
36 // ------------------------------------------- 36 // -------------------------------------------------------------------
37 37
38 #include <iomanip> 38 #include <iomanip>
39 39
40 #include "G4DNAMolecularReactionTable.hh" 40 #include "G4DNAMolecularReactionTable.hh"
41 #include "G4PhysicalConstants.hh" 41 #include "G4PhysicalConstants.hh"
42 #include "G4SystemOfUnits.hh" 42 #include "G4SystemOfUnits.hh"
43 #include "G4UIcommand.hh" 43 #include "G4UIcommand.hh"
44 #include "G4VDNAReactionModel.hh" 44 #include "G4VDNAReactionModel.hh"
45 #include "G4MoleculeHandleManager.hh" 45 #include "G4MoleculeHandleManager.hh"
46 #include "G4MoleculeTable.hh" 46 #include "G4MoleculeTable.hh"
47 #include "G4MolecularConfiguration.hh" 47 #include "G4MolecularConfiguration.hh"
48 #include "G4ReactionTableMessenger.hh" 48 #include "G4ReactionTableMessenger.hh"
49 #include "G4IosFlagsSaver.hh" 49 #include "G4IosFlagsSaver.hh"
50 #include "G4Exp.hh" 50 #include "G4Exp.hh"
51 51
52 using namespace std; 52 using namespace std;
53 53
54 G4DNAMolecularReactionTable* G4DNAMolecularRea << 54 G4DNAMolecularReactionTable* G4DNAMolecularReactionTable::fInstance(0);
55 55
56 G4DNAMolecularReactionData::G4DNAMolecularReac << 56 G4DNAMolecularReactionData::G4DNAMolecularReactionData() :
57 : fpReactant1(nullptr) << 57 fReactant1(),
58 , fpReactant2(nullptr) << 58 fReactant2(),
59 , fObservedReactionRate(0.) << 59 fObservedReactionRate(0.),
60 , fActivationRate(0.) << 60 fEffectiveReactionRadius(0.),
61 , fDiffusionRate(0.) << 61 fProducts(0)
62 , fOnsagerRadius(0.) << 62 {
63 , fReactionRadius(0.) << 63 fReactionID = 0;
64 , fEffectiveReactionRadius(0.) << 64 }
65 , fProbability(0.) << 65
66 , fType(0) << 66 //______________________________________________________________________________
67 , fReactionID(0) << 67
68 { << 68 G4DNAMolecularReactionData::
69 } << 69 G4DNAMolecularReactionData(G4double reactionRate,
70 << 70 G4MolecularConfiguration* reactant1,
71 G4DNAMolecularReactionData::G4DNAMolecularReac << 71 G4MolecularConfiguration* reactant2) :
72 << 72 fProducts(0)
73 << 73 {
74 : fpReactant1(pReactant1) << 74 fObservedReactionRate = reactionRate;
75 , fpReactant2(pReactant2) << 75 SetReactant1(reactant1);
76 , fObservedReactionRate(reactionRate) << 76 SetReactant2(reactant2);
77 , fActivationRate(0.) << 77
78 , fDiffusionRate(0.) << 78 G4double sumDiffCoeff(0.);
79 , fOnsagerRadius(0.) << 79
80 , fReactionRadius(0.) << 80 if (reactant1 == reactant2)
81 , fEffectiveReactionRadius(0.) << 81 {
82 , fProbability(0.) << 82 sumDiffCoeff = reactant1->GetDiffusionCoefficient();
83 , fType(0) << 83 fEffectiveReactionRadius = fObservedReactionRate
84 , fReactionID(0) << 84 / (4 * pi * sumDiffCoeff * Avogadro);
85 { << 85 }
86 ComputeEffectiveRadius(); << 86 else
87 } << 87 {
88 << 88 sumDiffCoeff = reactant1->GetDiffusionCoefficient()
89 G4DNAMolecularReactionData::G4DNAMolecularReac << 89 + reactant2->GetDiffusionCoefficient();
90 << 90 fEffectiveReactionRadius = fObservedReactionRate /
91 << 91 (4 * pi * sumDiffCoeff * Avogadro);
92 : fpReactant1(nullptr) << 92 }
93 , fpReactant2(nullptr) << 93 fReactionID = 0;
94 , fObservedReactionRate(reactionRate) << 94 }
95 , fActivationRate(0.) << 95
96 , fDiffusionRate(0.) << 96 //______________________________________________________________________________
97 , fOnsagerRadius(0.) << 97
98 , fReactionRadius(0.) << 98 G4DNAMolecularReactionData::
99 , fEffectiveReactionRadius(0.) << 99 G4DNAMolecularReactionData(G4double reactionRate,
100 , fProbability(0.) << 100 const G4String& reactant1,
101 , fType(0) << 101 const G4String& reactant2) :
102 , fReactionID(0) << 102 fProducts(0)
103 { << 103 {
104 SetReactant1(reactant1); << 104 fObservedReactionRate = reactionRate;
105 SetReactant2(reactant2); << 105 SetReactant1(reactant1);
106 ComputeEffectiveRadius(); << 106 SetReactant2(reactant2);
>> 107
>> 108 G4double sumDiffCoeff(0.);
>> 109
>> 110 if (fReactant1 == fReactant2)
>> 111 {
>> 112 sumDiffCoeff = fReactant1->GetDiffusionCoefficient();
>> 113 fEffectiveReactionRadius = fObservedReactionRate
>> 114 / (4 * pi * sumDiffCoeff * Avogadro);
>> 115 }
>> 116 else
>> 117 {
>> 118 sumDiffCoeff = fReactant1->GetDiffusionCoefficient()
>> 119 + fReactant2->GetDiffusionCoefficient();
>> 120 fEffectiveReactionRadius = fObservedReactionRate /
>> 121 (4 * pi * sumDiffCoeff * Avogadro);
>> 122 }
>> 123 fReactionID = 0;
107 } 124 }
108 125
109 G4DNAMolecularReactionData::~G4DNAMolecularRea 126 G4DNAMolecularReactionData::~G4DNAMolecularReactionData()
110 { 127 {
111 fProducts.clear(); << 128 if (fProducts)
112 } << 129 {
113 << 130 fProducts->clear();
114 void G4DNAMolecularReactionData::ComputeEffect << 131 delete fProducts;
115 { << 132 fProducts = 0;
116 G4double sumDiffCoeff = 0.; << 133 }
117 <<
118 if (fpReactant1 == fpReactant2) <<
119 { <<
120 sumDiffCoeff = fpReactant1->GetDiffusi <<
121 fEffectiveReactionRadius = fObservedRe <<
122 } <<
123 else <<
124 { <<
125 sumDiffCoeff = fpReactant1->GetDiffusi <<
126 + fpReactant2->GetDiffusi <<
127 fEffectiveReactionRadius = fObservedRe <<
128 } <<
129 <<
130 fReactionID = 0; <<
131 fReactionRadius = fEffectiveReactionRadius <<
132 fOnsagerRadius = (fpReactant1->GetCharge() <<
133 fProbability = 1; <<
134 <<
135 } 134 }
136 135
137 int G4DNAMolecularReactionData::GetReactionID( << 136 void G4DNAMolecularReactionData::SetReactant1(G4MolecularConfiguration* reactive)
138 { 137 {
139 return fReactionID; << 138 fReactant1 = reactive;
140 } 139 }
141 140
142 void G4DNAMolecularReactionData::SetReactionID <<
143 { <<
144 fReactionID = ID; <<
145 } <<
146 <<
147 void G4DNAMolecularReactionData::SetReactant1( <<
148 { <<
149 fpReactant1 = pReactive; <<
150 } <<
151 <<
152 void G4DNAMolecularReactionData::SetReactant2( <<
153 { <<
154 fpReactant2 = pReactive; <<
155 } <<
156 <<
157 void G4DNAMolecularReactionData::SetReactants( <<
158 <<
159 { <<
160 fpReactant1 = pReactant1; <<
161 fpReactant2 = pReactant2; <<
162 } <<
163 <<
164 void G4DNAMolecularReactionData::AddProduct(Re <<
165 { <<
166 fProducts.push_back(pMolecule); <<
167 } <<
168 141
169 G4int G4DNAMolecularReactionData::GetNbProduct << 142 void G4DNAMolecularReactionData::SetReactant2(G4MolecularConfiguration* reactive)
170 { 143 {
171 return (G4int)fProducts.size(); << 144 fReactant2 = reactive;
172 } 145 }
173 146
174 G4DNAMolecularReactionData::Reactant* G4DNAMol << 147 void G4DNAMolecularReactionData::SetReactants(G4MolecularConfiguration* reactant1,
>> 148 G4MolecularConfiguration* reactant2)
175 { 149 {
176 return fProducts[i]; << 150 fReactant1 = reactant1;
>> 151 fReactant2 = reactant2;
177 } 152 }
178 153
179 const G4DNAMolecularReactionData::ReactionProd << 154 void G4DNAMolecularReactionData::AddProduct(G4MolecularConfiguration* molecule)
180 { 155 {
181 return &fProducts; << 156 if (!fProducts) fProducts = new std::vector<G4MolecularConfiguration*>();
182 } << 157 fProducts->push_back(molecule);
183 <<
184 void G4DNAMolecularReactionData::RemoveProduct <<
185 { <<
186 fProducts.clear(); <<
187 } 158 }
188 159
189 void G4DNAMolecularReactionData::SetReactant1( 160 void G4DNAMolecularReactionData::SetReactant1(const G4String& reactive)
190 { 161 {
191 fpReactant1 = G4MoleculeTable::Instance()- << 162 fReactant1 = G4MoleculeTable::Instance()->GetConfiguration(reactive);
192 } 163 }
193 void G4DNAMolecularReactionData::SetReactant2( 164 void G4DNAMolecularReactionData::SetReactant2(const G4String& reactive)
194 { 165 {
195 fpReactant2 = G4MoleculeTable::Instance()- << 166 fReactant2 = G4MoleculeTable::Instance()->GetConfiguration(reactive);
196 } 167 }
197 void G4DNAMolecularReactionData::SetReactants( 168 void G4DNAMolecularReactionData::SetReactants(const G4String& reactant1,
198 << 169 const G4String& reactant2)
199 { 170 {
200 fpReactant1 = G4MoleculeTable::Instance()- << 171 fReactant1 = G4MoleculeTable::Instance()->GetConfiguration(reactant1);
201 fpReactant2 = G4MoleculeTable::Instance()- << 172 fReactant2 = G4MoleculeTable::Instance()->GetConfiguration(reactant2);
202 } 173 }
203 174
204 G4DNAMolecularReactionData::ReactantPair G4DNA << 175 void G4DNAMolecularReactionData::AddProduct(const G4String& molecule)
205 { 176 {
206 return std::make_pair(fpReactant1, fpReact << 177 if (!fProducts) fProducts = new std::vector<G4MolecularConfiguration*>();
>> 178 fProducts->push_back(G4MoleculeTable::Instance()->GetConfiguration(molecule));
207 } 179 }
208 180
209 G4DNAMolecularReactionData::Reactant* G4DNAMol <<
210 { <<
211 return fpReactant1; <<
212 } <<
213 181
214 G4DNAMolecularReactionData::Reactant* G4DNAMol << 182 double G4DNAMolecularReactionData::PolynomialParam(double temp_K, std::vector<double> P)
215 { << 183 {
216 return fpReactant2; << 184 double inv_temp = 1. / temp_K;
217 } <<
218 <<
219 void G4DNAMolecularReactionData::SetObservedRe <<
220 { <<
221 fObservedReactionRate = rate; <<
222 } <<
223 185
224 G4double G4DNAMolecularReactionData::GetObserv << 186 return pow(10,
225 { << 187 P[0] + P[1] * inv_temp + P[2] * pow(inv_temp, 2)
226 return fObservedReactionRate; << 188 + P[3] * pow(inv_temp, 3) + P[4] * pow(inv_temp, 4))
227 } << 189 * (1e-3 * CLHEP::m3 / (CLHEP::mole * CLHEP::s));
>> 190 }
>> 191
>> 192 double G4DNAMolecularReactionData::ArrehniusParam(double temp_K, std::vector<double> P)
>> 193 {
>> 194 return P[0]*G4Exp(P[1]/temp_K)*
>> 195 (1e-3 * CLHEP::m3 / (CLHEP::mole * CLHEP::s));
>> 196 }
>> 197
>> 198 double G4DNAMolecularReactionData::ScaledParameterization(double temp_K,
>> 199 double temp_init,
>> 200 double rateCste_init)
>> 201 {
>> 202 double D0 = G4MolecularConfiguration::DiffCoeffWater(temp_init);
>> 203 double Df = G4MolecularConfiguration::DiffCoeffWater(temp_K);
>> 204 return Df*rateCste_init/D0;
>> 205 }
228 206
229 G4double G4DNAMolecularReactionData::GetActiva << 207 //==============================================================================
230 { << 208 // REACTION TABLE
231 return fActivationRate; << 209 //==============================================================================
232 } <<
233 210
234 G4double G4DNAMolecularReactionData::GetDiffus << 211 G4DNAMolecularReactionTable* G4DNAMolecularReactionTable::GetReactionTable()
235 { 212 {
236 return fDiffusionRate; << 213 if (!fInstance)
>> 214 {
>> 215 fInstance = new G4DNAMolecularReactionTable();
>> 216 }
>> 217 return fInstance;
237 } 218 }
238 219
239 void G4DNAMolecularReactionData::SetReactionRa << 220 //_____________________________________________________________________________________
240 { <<
241 fReactionRadius = radius; <<
242 fEffectiveReactionRadius = -fOnsagerRadius <<
243 } <<
244 221
245 G4double G4DNAMolecularReactionData::GetReacti << 222 G4DNAMolecularReactionTable* G4DNAMolecularReactionTable::Instance()
246 { 223 {
247 return fReactionRadius; << 224 if (!fInstance)
>> 225 {
>> 226 fInstance = new G4DNAMolecularReactionTable();
>> 227 }
>> 228 return fInstance;
248 } 229 }
249 230
250 void G4DNAMolecularReactionData::SetEffectiveR << 231 //_____________________________________________________________________________________
251 { <<
252 fEffectiveReactionRadius = radius; <<
253 } <<
254 232
255 G4double G4DNAMolecularReactionData::GetEffect << 233 void G4DNAMolecularReactionTable::DeleteInstance()
256 { 234 {
257 return fEffectiveReactionRadius; << 235 // DEBUG
>> 236 // G4cout << "G4MolecularReactionTable::DeleteInstance" << G4endl;
>> 237 if (fInstance) delete fInstance;
>> 238 fInstance = 0;
258 } 239 }
259 240
260 G4double G4DNAMolecularReactionData::GetOnsage << 241 //_____________________________________________________________________________________
261 { <<
262 return fOnsagerRadius; <<
263 } <<
264 242
265 G4double G4DNAMolecularReactionData::GetProbab << 243 G4DNAMolecularReactionTable::G4DNAMolecularReactionTable() :
266 { << 244 G4ITReactionTable()
267 return fProbability; << 245 // , fMoleculeHandleManager(G4MoleculeHandleManager::Instance())
>> 246 {
>> 247 // G4cout << "G4DNAMolecularReactionTable::G4DNAMolecularReactionTable()" << G4endl;
>> 248 fVerbose = false;
>> 249 fpMessenger = new G4ReactionTableMessenger(this);
>> 250 return;
268 } 251 }
269 252
270 void G4DNAMolecularReactionData::SetProbabilit << 253 //_____________________________________________________________________________________
271 { <<
272 fProbability = prob; <<
273 } <<
274 254
275 void G4DNAMolecularReactionData::SetReactionTy << 255 G4DNAMolecularReactionTable::~G4DNAMolecularReactionTable()
276 { 256 {
277 G4double sumDiffCoeff = 0.; << 257 // DEBUG
>> 258 // G4cout << "G4MolecularReactionTable::~G4MolecularReactionTable" << G4endl;
278 259
279 if(type == 1) << 260 if(fpMessenger) delete fpMessenger;
280 { <<
281 261
282 sumDiffCoeff = fpReactant1->GetDiffusi << 262 ReactionDataMap::iterator it1 = fReactionData.begin();
283 fpReactant2->GetDiffusi << 263 std::map<G4MolecularConfiguration*,
>> 264 const G4DNAMolecularReactionData*>::iterator it2;
284 265
285 fReactionRadius = fpReactant1->GetVanD << 266 for(; it1 != fReactionData.end(); it1++)
286 fpReactant2->GetVanD << 267 {
>> 268 for(it2 = it1->second.begin(); it2 != it1->second.end(); it2++)
>> 269 {
>> 270 const G4DNAMolecularReactionData* reactionData = it2->second;
>> 271 if(reactionData)
>> 272 {
>> 273 G4MolecularConfiguration* reactant1 =
>> 274 reactionData->GetReactant1();
>> 275 G4MolecularConfiguration* reactant2 =
>> 276 reactionData->GetReactant2();
287 277
288 G4double Rs = 0.29 * nm; << 278 fReactionData[reactant1][reactant2] = 0;
289 if(fOnsagerRadius == 0) // Type II << 279 fReactionData[reactant2][reactant1] = 0;
290 { <<
291 fEffectiveReactionRadius = fReacti <<
292 fDiffusionRate = 4 * pi * sumDiffC <<
293 if (fpReactant1 == fpReactant2) fD <<
294 fActivationRate = fDiffusionRate * <<
295 fProbability = Rs / (Rs + (fDiffu <<
296 <<
297 }else{ // Type IV <<
298 fEffectiveReactionRadius = -fOnsag <<
299 fDiffusionRate = 4 * pi * sumDiffC <<
300 if (fpReactant1 == fpReactant2) fD <<
301 280
302 fActivationRate = fDiffusionRate * << 281 delete reactionData;
303 fProbability = Rs / (Rs + (fDiffus << 282 }
304 } <<
305 } 283 }
>> 284 }
306 285
307 fType = type; << 286 fReactionDataMV.clear();
308 } << 287 fReactionData.clear();
309 << 288 fReactantsMV.clear();
310 G4int G4DNAMolecularReactionData::GetReactionT <<
311 { <<
312 return fType; <<
313 } <<
314 <<
315 void G4DNAMolecularReactionData::AddProduct(co <<
316 { <<
317 fProducts.push_back(G4MoleculeTable::Insta <<
318 } <<
319 <<
320 double G4DNAMolecularReactionData::PolynomialP <<
321 { <<
322 double inv_temp = 1. / temp_K; <<
323 <<
324 return pow(10, <<
325 P[0] + P[1] * inv_temp + P[2] * <<
326 + P[3] * pow(inv_temp, 3) + P[4 <<
327 * (1e-3 * CLHEP::m3 / (CLHEP::mole * C <<
328 } <<
329 <<
330 double G4DNAMolecularReactionData::ArrehniusPa <<
331 { <<
332 return P[0] * G4Exp(P[1] / temp_K)* <<
333 (1e-3 * CLHEP::m3 / (CLHEP::mole * CLH <<
334 } <<
335 <<
336 double G4DNAMolecularReactionData::ScaledParam <<
337 <<
338 <<
339 { <<
340 double D0 = G4MolecularConfiguration::Diff <<
341 double Df = G4MolecularConfiguration::Diff <<
342 return Df * rateCste_init / D0; <<
343 } <<
344 <<
345 //============================================ <<
346 // REACTION TABLE <<
347 //============================================ <<
348 <<
349 G4DNAMolecularReactionTable* G4DNAMolecularRea <<
350 { <<
351 if (fpInstance == nullptr) <<
352 { <<
353 fpInstance = new G4DNAMolecularReactio <<
354 } <<
355 return fpInstance; <<
356 } 289 }
357 290
358 //____________________________________________ 291 //_____________________________________________________________________________________
359 292
360 G4DNAMolecularReactionTable* G4DNAMolecularRea << 293 void G4DNAMolecularReactionTable::SetReaction(G4DNAMolecularReactionData* reactionData)
361 { 294 {
362 if (fpInstance == nullptr) << 295 G4MolecularConfiguration* reactant1 = reactionData->GetReactant1();
363 { << 296 G4MolecularConfiguration* reactant2 = reactionData->GetReactant2();
364 fpInstance = new G4DNAMolecularReactio <<
365 } <<
366 return fpInstance; <<
367 } <<
368 297
369 //____________________________________________ << 298 fReactionData[reactant1][reactant2] = reactionData;
>> 299 fReactantsMV[reactant1].push_back(reactant2);
>> 300 fReactionDataMV[reactant1].push_back(reactionData);
>> 301
>> 302 if (reactant1 != reactant2)
>> 303 {
>> 304 fReactionData[reactant2][reactant1] = reactionData;
>> 305 fReactantsMV[reactant2].push_back(reactant1);
>> 306 fReactionDataMV[reactant2].push_back(reactionData);
>> 307 }
370 308
371 void G4DNAMolecularReactionTable::DeleteInstan << 309 fVectorOfReactionData.push_back(reactionData);
372 { <<
373 <<
374 <<
375 delete fpInstance; <<
376 <<
377 fpInstance = nullptr; <<
378 } 310 }
379 311
380 //____________________________________________ 312 //_____________________________________________________________________________________
381 313
382 G4DNAMolecularReactionTable::G4DNAMolecularRea << 314 void G4DNAMolecularReactionTable::SetReaction(G4double reactionRate,
383 : << 315 G4MolecularConfiguration* reactant1,
384 fpMessenger(new G4ReactionTableMessenger( << 316 G4MolecularConfiguration* reactant2)
385 { 317 {
>> 318 G4DNAMolecularReactionData* reactionData = new G4DNAMolecularReactionData(
>> 319 reactionRate, reactant1, reactant2);
>> 320 SetReaction(reactionData);
386 } 321 }
387 322
388 G4DNAMolecularReactionTable::~G4DNAMolecularRe << 323 //_____________________________________________________________________________________
389 324
390 void G4DNAMolecularReactionTable::SetReaction( << 325 void G4DNAMolecularReactionTable::PrintTable(G4VDNAReactionModel* pReactionModel)
391 { 326 {
392 const auto pReactant1 = pReactionData->Get << 327 // Print Reactions and Interaction radius for jump step = 3ps
393 const auto pReactant2 = pReactionData->Get <<
394 <<
395 fReactionData[pReactant1][pReactant2] = pR <<
396 fReactantsMV[pReactant1].push_back(pReacta <<
397 fReactionDataMV[pReactant1].push_back(pRea <<
398 328
399 if (pReactant1 != pReactant2) << 329 G4IosFlagsSaver iosfs(G4cout);
400 { <<
401 fReactionData[pReactant2][pReactant1] <<
402 fReactantsMV[pReactant2].push_back(pRe <<
403 fReactionDataMV[pReactant2].push_back( <<
404 } <<
405 330
406 fVectorOfReactionData.emplace_back(pReacti << 331 if (pReactionModel)
407 pReactionData->SetReactionID((G4int)fVecto << 332 {
408 } << 333 if (!(pReactionModel->GetReactionTable())) pReactionModel->SetReactionTable(
>> 334 this);
>> 335 }
409 336
410 //____________________________________________ << 337 ReactivesMV::iterator itReactives;
411 338
412 void G4DNAMolecularReactionTable::SetReaction( << 339 map<G4MolecularConfiguration*, map<G4MolecularConfiguration*, G4bool> > alreadyPrint;
413 <<
414 <<
415 { <<
416 auto reactionData = new G4DNAMolecularReac <<
417 SetReaction(reactionData); <<
418 } <<
419 340
420 //____________________________________________ << 341 G4cout << "Number of chemical species involved in reactions = "
>> 342 << fReactantsMV.size() << G4endl;
421 343
422 void G4DNAMolecularReactionTable::PrintTable(G << 344 G4int nbPrintable = fReactantsMV.size() * fReactantsMV.size();
423 { <<
424 // Print Reactions and Interaction radius <<
425 345
426 G4IosFlagsSaver iosfs(G4cout); << 346 G4String *outputReaction = new G4String[nbPrintable];
>> 347 G4String *outputReactionRate = new G4String[nbPrintable];
>> 348 G4String *outputRange = new G4String[nbPrintable];
>> 349 G4int n = 0;
427 350
428 if ((pReactionModel != nullptr) && ((pReact << 351 for (itReactives = fReactantsMV.begin(); itReactives != fReactantsMV.end();
429 { << 352 itReactives++)
430 pReactionModel->SetReactionTable(this); << 353 {
431 } << 354 G4MolecularConfiguration* moleculeA = (G4MolecularConfiguration*) itReactives->first;
>> 355 const vector<G4MolecularConfiguration*>* reactivesVector = CanReactWith(moleculeA);
432 356
433 ReactivesMV::iterator itReactives; << 357 if (pReactionModel) pReactionModel->InitialiseToPrint(moleculeA);
434 358
435 std::map<Reactant*, std::map<Reactant*, G4 << 359 G4int nbReactants = fReactantsMV[itReactives->first].size();
436 360
437 G4cout << "Number of chemical species invo << 361 for (G4int iReact = 0; iReact < nbReactants; iReact++)
438 << fReactantsMV.size() << G4endl; << 362 {
439 363
440 std::size_t nbPrintable = fReactantsMV.siz << 364 G4MolecularConfiguration* moleculeB = (G4MolecularConfiguration*) (*reactivesVector)[iReact];
441 365
442 auto outputReaction = new G4String[nbPrin << 366 const G4DNAMolecularReactionData* reactionData =
443 auto outputReactionRate = new G4String[nb << 367 fReactionData[moleculeA][moleculeB];
444 auto outputRange = new G4String[nbPrintab <<
445 G4int n = 0; <<
446 368
447 for (itReactives = fReactantsMV.begin(); i << 369 //-----------------------------------------------------------
448 ++itReactives) << 370 // Name of the reaction
449 { << 371 if (!alreadyPrint[moleculeA][moleculeB])
450 auto moleculeA = (Reactant*)itReactiv << 372 {
451 const vector<Reactant*>* reactivesVect << 373 outputReaction[n] = moleculeA->GetName() + " + " + moleculeB->GetName();
452 374
453 if (pReactionModel != nullptr) pReacti << 375 G4int nbProducts = reactionData->GetNbProducts();
454 376
455 auto nbReactants = (G4int)fReactantsM << 377 if (nbProducts)
>> 378 {
>> 379 outputReaction[n] += " -> " + reactionData->GetProduct(0)->GetName();
456 380
457 for (G4int iReact = 0; iReact < nbReac << 381 for (G4int j = 1; j < nbProducts; j++)
>> 382 {
>> 383 outputReaction[n] += " + " + reactionData->GetProduct(j)->GetName();
>> 384 }
>> 385 }
>> 386 else
458 { 387 {
459 auto moleculeB = (Reactant*)(*reac << 388 outputReaction[n] += " -> No product";
>> 389 }
460 390
461 Data* reactionData = fReactionData << 391 //-----------------------------------------------------------
>> 392 // Interaction Rate
>> 393 outputReactionRate[n] = G4UIcommand::ConvertToString(
>> 394 reactionData->GetObservedReactionRateConstant() / (1e-3 * m3 / (mole * s)));
>> 395
>> 396 //-----------------------------------------------------------
>> 397 // Calculation of the Interaction Range
>> 398 G4double interactionRange = -1;
>> 399 if (pReactionModel) interactionRange =
>> 400 pReactionModel->GetReactionRadius(iReact);
462 401
463 //-------------------------------- << 402 if (interactionRange != -1)
464 // Name of the reaction << 403 {
465 if (!alreadyPrint[moleculeA][molec << 404 outputRange[n] = G4UIcommand::ConvertToString(
466 { << 405 interactionRange / nanometer);
467 outputReaction[n] = moleculeA- << 406 }
468 << 407 else
469 G4int nbProducts = reactionDat << 408 {
470 << 409 outputRange[n] = "";
471 if (nbProducts != 0) <<
472 { <<
473 outputReaction[n] += " -> <<
474 <<
475 for (G4int j = 1; j < nbPr <<
476 { <<
477 outputReaction[n] += " <<
478 } <<
479 } <<
480 else <<
481 { <<
482 outputReaction[n] += " -> <<
483 } <<
484 <<
485 //---------------------------- <<
486 // Interaction Rate <<
487 outputReactionRate[n] = G4UIco <<
488 reactionData->GetObservedR <<
489 <<
490 //---------------------------- <<
491 // Calculation of the Interact <<
492 G4double interactionRange = -1 <<
493 if (pReactionModel != nullptr) <<
494 pReactionModel->GetReactio <<
495 <<
496 if (interactionRange != -1) <<
497 { <<
498 outputRange[n] = G4UIcomma <<
499 interactionRange / nan <<
500 } <<
501 else <<
502 { <<
503 outputRange[n] = ""; <<
504 } <<
505 <<
506 alreadyPrint[moleculeB][molecu <<
507 n++; <<
508 } <<
509 } 410 }
>> 411
>> 412 alreadyPrint[moleculeB][moleculeA] = TRUE;
>> 413 n++;
>> 414 }
510 } 415 }
511 // G4cout<<"Number of possible reactions: << 416 }
>> 417 // G4cout<<"Number of possible reactions: "<< n << G4endl;
512 418
513 ////////////////////////////////////////// << 419 ////////////////////////////////////////////////////////////////////
514 // Tableau dynamique en fonction du nombre << 420 // Tableau dynamique en fonction du nombre de caractere maximal dans
515 // chaque colonne << 421 // chaque colonne
516 ////////////////////////////////////////// << 422 ////////////////////////////////////////////////////////////////////
517 423
518 G4int maxlengthOutputReaction = -1; << 424 G4int maxlengthOutputReaction = -1;
519 G4int maxlengthOutputReactionRate = -1; << 425 G4int maxlengthOutputReactionRate = -1;
520 426
521 for (G4int i = 0; i < n; ++i) << 427 for (G4int i = 0; i < n; i++)
>> 428 {
>> 429 if (maxlengthOutputReaction < (G4int) outputReaction[i].length())
522 { 430 {
523 if (maxlengthOutputReaction < (G4int)o << 431 maxlengthOutputReaction = outputReaction[i].length();
524 { << 432 }
525 maxlengthOutputReaction = (G4int)o << 433 if (maxlengthOutputReactionRate < (G4int) outputReactionRate[i].length())
526 } << 434 {
527 if (maxlengthOutputReactionRate < (G4i << 435 maxlengthOutputReactionRate = outputReactionRate[i].length();
528 { <<
529 maxlengthOutputReactionRate = (G4i <<
530 } <<
531 } 436 }
>> 437 }
>> 438
>> 439 maxlengthOutputReaction += 2;
>> 440 maxlengthOutputReactionRate += 2;
>> 441
>> 442 if (maxlengthOutputReaction < 10) maxlengthOutputReaction = 10;
>> 443 if (maxlengthOutputReactionRate < 30) maxlengthOutputReactionRate = 30;
532 444
533 maxlengthOutputReaction += 2; << 445 G4String* title;
534 maxlengthOutputReactionRate += 2; <<
535 446
536 if (maxlengthOutputReaction < 10) maxlengt << 447 if (pReactionModel) title = new G4String[3];
537 if (maxlengthOutputReactionRate < 30) maxl << 448 else title = new G4String[2];
538 449
539 G4String* title; << 450 title[0] = "Reaction";
>> 451 title[1] = "Reaction Rate [dm3/(mol*s)]";
540 452
541 if (pReactionModel != nullptr) title = new << 453 if (pReactionModel) title[2] =
542 else title = new G4String[2]; << 454 "Interaction Range for chosen reaction model [nm]";
543 455
544 title[0] = "Reaction"; << 456 G4cout << setfill(' ') << setw(maxlengthOutputReaction) << left << title[0]
545 title[1] = "Reaction Rate [dm3/(mol*s)]"; << 457 << setw(maxlengthOutputReactionRate) << left << title[1];
546 458
547 if (pReactionModel != nullptr) title[2] = << 459 if (pReactionModel) G4cout << setw(2) << left << title[2];
548 "Interaction Range for chosen reaction <<
549 460
550 G4cout << setfill(' ') << setw(maxlengthOu << 461 G4cout << G4endl;
551 << setw(maxlengthOutputReactionRate) < <<
552 462
553 if (pReactionModel != nullptr) G4cout << s << 463 G4cout.fill('-');
>> 464 if (pReactionModel) G4cout.width(
>> 465 maxlengthOutputReaction + 2 + maxlengthOutputReactionRate + 2
>> 466 + (G4int) title[2].length());
>> 467 else G4cout.width(maxlengthOutputReaction + 2 + maxlengthOutputReactionRate);
>> 468 G4cout << "-" << G4endl;
>> 469 G4cout.fill(' ');
>> 470
>> 471 for (G4int i = 0; i < n; i++)
>> 472 {
>> 473 G4cout << setw(maxlengthOutputReaction) << left << outputReaction[i]
>> 474 << setw(maxlengthOutputReactionRate) << left
>> 475 << outputReactionRate[i];
>> 476
>> 477 if (pReactionModel) G4cout << setw(2) << left << outputRange[i];
554 478
555 G4cout << G4endl; 479 G4cout << G4endl;
556 480
557 G4cout.fill('-'); 481 G4cout.fill('-');
558 if (pReactionModel != nullptr) G4cout.widt << 482 if (pReactionModel) G4cout.width(
559 maxlengthOutputReaction + 2 + maxlengt 483 maxlengthOutputReaction + 2 + maxlengthOutputReactionRate + 2
560 + (G4int)title[2].length()); << 484 + (G4int) title[2].length());
561 else G4cout.width(maxlengthOutputReaction << 485 else G4cout.width(
>> 486 maxlengthOutputReaction + 2 + maxlengthOutputReactionRate);
562 G4cout << "-" << G4endl; 487 G4cout << "-" << G4endl;
563 G4cout.fill(' '); 488 G4cout.fill(' ');
>> 489 }
564 490
565 for (G4int i = 0; i < n; i++) << 491 delete[] title;
566 { << 492 delete[] outputReaction;
567 G4cout << setw(maxlengthOutputReaction << 493 delete[] outputReactionRate;
568 << setw(maxlengthOutputReaction << 494 delete[] outputRange;
569 << outputReactionRate[i]; <<
570 <<
571 if (pReactionModel != nullptr) G4cout <<
572 <<
573 G4cout << G4endl; <<
574 <<
575 G4cout.fill('-'); <<
576 if (pReactionModel != nullptr) G4cout. <<
577 maxlengthOutputReaction + 2 + maxl <<
578 + (G4int)title[2].length()); <<
579 else G4cout.width( <<
580 maxlengthOutputReaction + 2 + maxl <<
581 G4cout << "-" << G4endl; <<
582 G4cout.fill(' '); <<
583 } <<
584 <<
585 delete[] title; <<
586 delete[] outputReaction; <<
587 delete[] outputReactionRate; <<
588 delete[] outputRange; <<
589 } 495 }
590 496
591 //____________________________________________ 497 //______________________________________________________________________________
592 // Get/Set methods 498 // Get/Set methods
593 499
594 G4VDNAMolecularGeometry* G4DNAMolecularReactio << 500 const G4DNAMolecularReactionData*
595 { << 501 G4DNAMolecularReactionTable::GetReactionData(G4MolecularConfiguration* reactant1,
596 return fGeometry; << 502 G4MolecularConfiguration* reactant2) const
597 } << 503 {
598 << 504 if (fReactionData.empty())
599 G4DNAMolecularReactionTable::Data* << 505 {
600 G4DNAMolecularReactionTable::GetReactionData(R << 506 G4String errMsg = "No reaction table was implemented";
601 R << 507 G4Exception("G4MolecularInteractionTable::GetReactionData", "",
602 { << 508 FatalErrorInArgument, errMsg);
603 if (fReactionData.empty()) << 509 return 0;
604 { << 510 }
605 G4String errMsg = "No reaction table w << 511
606 G4Exception("G4MolecularInteractionTab << 512 ReactionDataMap::const_iterator it1 = fReactionData.find(reactant1);
607 FatalErrorInArgument, errM << 513
608 } << 514 if (it1 == fReactionData.end())
609 << 515 {
610 auto it1 = fReactionData.find(pReactant1); << 516 G4String errMsg =
611 << 517 "No reaction table was implemented for this molecule Definition : " + reactant1
612 if (it1 == fReactionData.end()) <<
613 { <<
614 G4String errMsg = <<
615 "No reaction table was implemented <<
616 ->GetName(); 518 ->GetName();
617 G4Exception("G4MolecularInteractionTab << 519 // G4cout << "--- G4MolecularInteractionTable::GetReactionData ---" << G4endl;
618 FatalErrorInArgument, errM << 520 // G4cout << errMsg << G4endl;
619 // Though the above is Fatal and will << 521 G4Exception("G4MolecularInteractionTable::GetReactionData", "",
620 return nullptr; << 522 FatalErrorInArgument, errMsg);
621 } << 523 // return 0;
>> 524 }
622 525
623 auto it2 = it1->second.find(pReactant2); << 526 std::map<G4MolecularConfiguration*,
>> 527 const G4DNAMolecularReactionData*>::const_iterator it2 =
>> 528 it1->second.find(reactant2);
624 529
625 if (it2 == it1->second.end()) << 530 if (it2 == it1->second.end())
626 { << 531 {
627 G4cout << "Name : " << pReactant2->Get << 532 G4cout << "Name : " << reactant2->GetName() << G4endl;
628 G4String errMsg = "No reaction table w << 533 G4String errMsg = "No reaction table was implemented for this molecule : "
629 + pReactant2->GetName(); << 534 + reactant2 -> GetName();
630 G4Exception("G4MolecularInteractionTab << 535 G4Exception("G4MolecularInteractionTable::GetReactionData","",FatalErrorInArgument, errMsg);
631 } << 536 }
632 537
633 return (it2->second); << 538 return (it2->second);
634 } <<
635 <<
636 const G4DNAMolecularReactionTable::ReactionDat <<
637 { <<
638 return fReactionData; <<
639 } <<
640 <<
641 G4DNAMolecularReactionTable::DataList G4DNAMol <<
642 { <<
643 DataList dataList; <<
644 <<
645 for (const auto& pData : fVectorOfReaction <<
646 { <<
647 dataList.emplace_back(pData.get()); <<
648 } <<
649 <<
650 return dataList; <<
651 } 539 }
652 540
653 //____________________________________________ 541 //______________________________________________________________________________
654 542
655 const G4DNAMolecularReactionTable::ReactantLis << 543 const std::vector<G4MolecularConfiguration*>*
656 G4DNAMolecularReactionTable::CanReactWith(Reac << 544 G4DNAMolecularReactionTable::CanReactWith(G4MolecularConfiguration * aMolecule) const
657 { 545 {
658 if (fReactantsMV.empty()) << 546 if (fReactantsMV.empty())
659 { << 547 {
660 G4String errMsg = "No reaction table w << 548 G4String errMsg = "No reaction table was implemented";
661 G4Exception("G4MolecularInteractionTab << 549 G4Exception("G4MolecularInteractionTable::CanReactWith", "",
662 FatalErrorInArgument, errM << 550 FatalErrorInArgument, errMsg);
663 return nullptr; << 551 return 0;
664 } << 552 }
665 553
666 auto itReactivesMap = fReactantsMV.find(pM << 554 ReactivesMV::const_iterator itReactivesMap = fReactantsMV.find(aMolecule);
667 555
668 if (itReactivesMap == fReactantsMV.end()) << 556 if (itReactivesMap == fReactantsMV.end())
669 { << 557 {
670 #ifdef G4VERBOSE 558 #ifdef G4VERBOSE
671 if (fVerbose) <<
672 { <<
673 G4String errMsg = "No reaction tab <<
674 + pMolecule->GetName(); <<
675 // G4Exception("G4Molecul <<
676 G4cout << "--- G4MolecularInteract <<
677 G4cout << errMsg << G4endl; <<
678 } <<
679 #endif <<
680 return nullptr; <<
681 } <<
682 <<
683 if (fVerbose) 559 if (fVerbose)
684 { 560 {
685 G4cout << " G4MolecularInteractionTabl << 561 G4String errMsg = "No reaction table was implemented for this molecule : "
686 G4cout << "You are checking reactants << 562 + aMolecule->GetName();
687 G4cout << " the number of reactants is << 563 // G4Exception("G4MolecularInteractionTable::CanReactWith","",FatalErrorInArgument, errMsg);
688 << 564 G4cout << "--- G4MolecularInteractionTable::GetReactionData ---" << G4endl;
689 auto itProductsVector = itReactivesMap << 565 G4cout << errMsg << G4endl;
690 << 566 }
691 for (; itProductsVector != itReactives << 567 #endif
692 { << 568 return 0;
693 G4cout << (*itProductsVector)->Get << 569 }
694 } << 570 else
>> 571 {
>> 572 if(fVerbose)
>> 573 {
>> 574 G4cout<< " G4MolecularInteractionTable::CanReactWith :"<<G4endl;
>> 575 G4cout<<"You are checking reactants for : " << aMolecule->GetName()<<G4endl;
>> 576 G4cout<<" the number of reactants is : " << itReactivesMap->second.size()<<G4endl;
>> 577
>> 578 std::vector<G4MolecularConfiguration*>::const_iterator itProductsVector =
>> 579 itReactivesMap->second.begin();
>> 580
>> 581 for(; itProductsVector != itReactivesMap->second.end(); itProductsVector++)
>> 582 {
>> 583 G4cout<<(*itProductsVector)->GetName()<<G4endl;
>> 584 }
695 } 585 }
696 return &(itReactivesMap->second); 586 return &(itReactivesMap->second);
>> 587 }
>> 588 return 0;
697 } 589 }
698 590
699 //____________________________________________ 591 //______________________________________________________________________________
700 592
701 const G4DNAMolecularReactionTable::SpecificDat << 593 const std::map<G4MolecularConfiguration*, const G4DNAMolecularReactionData*>*
702 G4DNAMolecularReactionTable::GetReativesNData( << 594 G4DNAMolecularReactionTable::GetReativesNData(G4MolecularConfiguration* molecule) const
703 { 595 {
704 if (fReactionData.empty()) << 596 if (fReactionData.empty())
705 { << 597 {
706 G4String errMsg = "No reaction table w << 598 G4String errMsg = "No reaction table was implemented";
707 G4Exception("G4MolecularInteractionTab << 599 G4Exception("G4MolecularInteractionTable::CanInteractWith", "",
708 FatalErrorInArgument, errM << 600 FatalErrorInArgument, errMsg);
709 } << 601 return 0;
710 << 602 }
711 auto itReactivesMap = fReactionData.find(m << 603
712 << 604 ReactionDataMap::const_iterator itReactivesMap = fReactionData.find(molecule);
713 if (itReactivesMap == fReactionData.end()) << 605
714 { << 606 if (itReactivesMap == fReactionData.end())
715 return nullptr; << 607 {
716 } << 608 return 0;
717 << 609 // G4cout << "Nom : " << molecule->GetName() << G4endl;
718 if (fVerbose) << 610 // G4String errMsg = "No reaction table was implemented for this molecule Definition : "
719 { << 611 // + molecule -> GetName();
720 G4cout << " G4MolecularInteractionTabl << 612 // G4Exception("G4MolecularInteractionTable::CanReactWith","",FatalErrorInArgument, errMsg);
721 G4cout << "You are checking reactants << 613 }
722 G4cout << " the number of reactants is << 614 else
723 << 615 {
724 auto itProductsVector = itReactivesMap << 616 if(fVerbose)
725 << 617 {
726 for (; itProductsVector != itReactives << 618 G4cout<< " G4MolecularInteractionTable::CanReactWith :"<<G4endl;
727 { << 619 G4cout<<"You are checking reactants for : " << molecule->GetName()<<G4endl;
728 G4cout << itProductsVector->first- << 620 G4cout<<" the number of reactants is : " << itReactivesMap->second.size()<<G4endl;
729 } << 621
>> 622 std::map<G4MolecularConfiguration*,
>> 623 const G4DNAMolecularReactionData*>::const_iterator itProductsVector =
>> 624 itReactivesMap->second.begin();
>> 625
>> 626 for(; itProductsVector != itReactivesMap->second.end(); itProductsVector++)
>> 627 {
>> 628 G4cout<<itProductsVector->first->GetName()<<G4endl;
>> 629 }
730 } 630 }
731 return &(itReactivesMap->second); 631 return &(itReactivesMap->second);
>> 632 }
>> 633
>> 634 return 0;
732 } 635 }
733 636
734 //____________________________________________ 637 //______________________________________________________________________________
735 638
736 const G4DNAMolecularReactionTable::DataList* << 639 const std::vector<const G4DNAMolecularReactionData*>*
737 G4DNAMolecularReactionTable::GetReactionData(c << 640 G4DNAMolecularReactionTable::GetReactionData(G4MolecularConfiguration* molecule) const
738 { 641 {
739 if (fReactionDataMV.empty()) << 642 if (fReactionDataMV.empty())
740 { << 643 {
741 G4String errMsg = "No reaction table w << 644 G4String errMsg = "No reaction table was implemented";
742 G4Exception("G4MolecularInteractionTab << 645 G4Exception("G4MolecularInteractionTable::CanInteractWith", "",
743 FatalErrorInArgument, errM << 646 FatalErrorInArgument, errMsg);
744 } << 647 return 0;
745 auto it = fReactionDataMV.find(molecule); << 648 }
>> 649 ReactionDataMV::const_iterator it = fReactionDataMV.find(molecule);
746 650
747 if (it == fReactionDataMV.end()) << 651 if (it == fReactionDataMV.end())
748 { << 652 {
749 G4String errMsg = "No reaction table w << 653 G4cout << "Nom : " << molecule->GetName() << G4endl;
750 + molecule->GetName(); << 654 G4String errMsg = "No reaction table was implemented for this molecule Definition : "
751 G4Exception("G4MolecularInteractionTab << 655 + molecule -> GetName();
752 // Though the above is Fatal and will << 656 G4Exception("G4MolecularInteractionTable::GetReactionData","",FatalErrorInArgument, errMsg);
753 return nullptr; << 657 return 0; // coverity
754 } << 658 }
755 659
756 return &(it->second); << 660 return &(it->second);
757 } 661 }
758 662
759 //____________________________________________ 663 //______________________________________________________________________________
760 664
761 G4DNAMolecularReactionTable::Data* G4DNAMolecu << 665 const G4DNAMolecularReactionData*
762 << 666 G4DNAMolecularReactionTable::GetReactionData(const G4String& mol1,
>> 667 const G4String& mol2) const
763 { 668 {
764 const auto pConf1 = G4MoleculeTable::GetMo << 669 G4MolecularConfiguration* conf1 = G4MoleculeTable::GetMoleculeTable()
765 const auto pConf2 = G4MoleculeTable::GetMo << 670 ->GetConfiguration(mol1);
766 return GetReactionData(pConf1, pConf2); << 671 G4MolecularConfiguration* conf2 = G4MoleculeTable::GetMoleculeTable()
>> 672 ->GetConfiguration(mol2);
>> 673
>> 674 return GetReactionData(conf1, conf2);
767 } 675 }
768 676
769 //____________________________________________ 677 //______________________________________________________________________________
770 678
771 void 679 void
772 G4DNAMolecularReactionData::SetPolynomialParam << 680 G4DNAMolecularReactionData::
>> 681 SetPolynomialParameterization(const std::vector<double>& P)
773 { 682 {
774 fRateParam = std::bind(PolynomialParam, st << 683 fRateParam = std::bind(PolynomialParam, std::placeholders::_1, P);
775 } 684 }
776 685
777 //____________________________________________ 686 //______________________________________________________________________________
778 687
779 void G4DNAMolecularReactionData::SetArrehniusP 688 void G4DNAMolecularReactionData::SetArrehniusParameterization(double A0,
780 689 double E_R)
781 { 690 {
782 std::vector<double> P = { A0, E_R }; << 691 std::vector<double> P = {A0, E_R};
783 fRateParam = std::bind(ArrehniusParam, std << 692
>> 693 G4cout << "ici = " << P[0] << G4endl;
>> 694 G4cout << "A0 = " << A0 << G4endl;
>> 695
>> 696 fRateParam = std::bind(ArrehniusParam, std::placeholders::_1, P);
784 } 697 }
785 698
786 //____________________________________________ 699 //______________________________________________________________________________
787 700
788 void G4DNAMolecularReactionData::SetScaledPara 701 void G4DNAMolecularReactionData::SetScaledParameterization(double temperature_K,
789 702 double rateCste)
790 { 703 {
791 fRateParam = std::bind(ScaledParameterizat << 704 fRateParam = std::bind(ScaledParameterization,
792 std::placeholders:: << 705 std::placeholders::_1,
793 temperature_K, << 706 temperature_K,
794 rateCste); << 707 rateCste);
795 } 708 }
796 709
797 //____________________________________________ 710 //______________________________________________________________________________
798 711
799 void G4DNAMolecularReactionTable::ScaleReactio 712 void G4DNAMolecularReactionTable::ScaleReactionRateForNewTemperature(double temp_K)
800 { 713 {
801 for (const auto& pData : fVectorOfReaction << 714 size_t end = fVectorOfReactionData.size();
802 { << 715
803 const_cast<G4DNAMolecularReactionData* << 716 for(size_t i = 0 ; i < end ; ++i)
804 } << 717 {
>> 718 ((G4DNAMolecularReactionData*) fVectorOfReactionData[i])->
>> 719 ScaleForNewTemperature(temp_K);
>> 720 }
805 } 721 }
806 722
807 //____________________________________________ 723 //______________________________________________________________________________
808 724
809 void G4DNAMolecularReactionData::ScaleForNewTe 725 void G4DNAMolecularReactionData::ScaleForNewTemperature(double temp_K)
810 { 726 {
811 if (fRateParam) << 727 if(fRateParam)
812 { << 728 {
813 SetObservedReactionRateConstant(fRateP << 729 SetObservedReactionRateConstant(fRateParam(temp_K));
814 } << 730
>> 731 // G4cout <<"PROD RATE = " << fProductionRate << G4endl;
>> 732 //
>> 733 // if(fProductionRate != DBL_MAX && fProductionRate !=0)
>> 734 // {
>> 735 // SetPartiallyDiffusionControlledReactionByActivation(fObservedReactionRate,
>> 736 // fProductionRate);
>> 737 // }
>> 738 }
815 } 739 }
816 740
817 //____________________________________________ 741 //______________________________________________________________________________
818 742
819 G4DNAMolecularReactionTable::Data* << 743 const G4DNAMolecularReactionData*
820 G4DNAMolecularReactionTable::GetReaction(int r 744 G4DNAMolecularReactionTable::GetReaction(int reactionID) const
821 { 745 {
822 for (auto& pData : fVectorOfReactionData) << 746 for(size_t i = 0 ; i < fVectorOfReactionData.size() ; ++i)
823 { << 747 {
824 if (pData->GetReactionID() == reaction << 748 if(fVectorOfReactionData[i]->GetReactionID() == reactionID)
825 { << 749 return fVectorOfReactionData[i];
826 return pData.get(); << 750 }
827 } << 751 return nullptr;
828 } <<
829 return nullptr; <<
830 } <<
831 <<
832 size_t G4DNAMolecularReactionTable::GetNReacti <<
833 { <<
834 return fVectorOfReactionData.size(); <<
835 } <<
836 <<
837 void G4DNAMolecularReactionTable::Reset() <<
838 { <<
839 fReactionData.clear(); <<
840 fReactantsMV.clear(); <<
841 fReactionDataMV.clear(); <<
842 fVectorOfReactionData.clear(); <<
843 } 752 }
844 753