diff --git a/source/materials/History b/source/materials/History index 7f90b0fb9f..ed14b01bf1 100644 --- a/source/materials/History +++ b/source/materials/History @@ -16,6 +16,19 @@ committal in the CVS repository ! * Reverse chronological order (last date on top), please * ---------------------------------------------------------- +## 2022-04-14 Vladimir Ivanchenko (materials-V10-07-28) +## 2022-04-12 Vladimir Ivanchenko +- G4DensityEffectData, G4IonisParamMat - fixed density effect correction + for liquid hydrogen (problem #2346), code clean-up +- G4NistMaterialBuilder - fixed mean ionisation potentilal of carbon + to the current NIST value from 81 eV to 78 eV (problem #2474) +- G4Material - reorganized AddElement(..) and AddMaterial(..) methods, + temporary vector are used that allows addition of the same elements + during initialisation, when addition completed the vector of elements + and arrays of fractions are created, filled and temporary vectors + are deleted, as a result the code become more transparent and + fixing problem #2486 + 04-02-22 V.Ivanchenko (materials-V10-07-27) - G4SandiaTable - fixed compilation warning on unused variable diff --git a/source/materials/include/G4DensityEffectData.hh b/source/materials/include/G4DensityEffectData.hh index 037447d0f1..0dd53648d7 100644 --- a/source/materials/include/G4DensityEffectData.hh +++ b/source/materials/include/G4DensityEffectData.hh @@ -61,7 +61,7 @@ public: ~G4DensityEffectData(); // return index by Z, -1 if material is not in the table - G4int GetElementIndex(G4int Z, G4State mState) const; + G4int GetElementIndex(G4int Z, G4State st = kStateUndefined) const; // return index by material name, -1 if material is not in the table G4int GetIndex(const G4String& matName) const; diff --git a/source/materials/include/G4Material.hh b/source/materials/include/G4Material.hh index a9a1b45ca5..04ac42a7b2 100644 --- a/source/materials/include/G4Material.hh +++ b/source/materials/include/G4Material.hh @@ -351,6 +351,11 @@ private: G4int fIdxComponent; // Index of a new component G4bool fMassFraction; // Flag of the method to add components + // For composites built + std::vector* fAtoms; + std::vector* fElmFrac; + std::vector* fElm; + // For composites built via AddMaterial() std::map fMatComponents; diff --git a/source/materials/src/G4DensityEffectData.cc b/source/materials/src/G4DensityEffectData.cc index 39493c75cd..ad4fe63d76 100644 --- a/source/materials/src/G4DensityEffectData.cc +++ b/source/materials/src/G4DensityEffectData.cc @@ -1294,13 +1294,9 @@ void G4DensityEffectData::Initialize() AddMaterial(M277,"G4_GRAPHITE_POROUS"); } -G4int G4DensityEffectData::GetElementIndex(G4int Z, G4State st) const +G4int G4DensityEffectData::GetElementIndex(G4int Z, G4State) const { - G4int idx = -1; - if(Z > 0 && Z < NDENSELEM) { - if(st == state[Z] || st == kStateUndefined) { idx = indexZ[Z]; } - } - return idx; + return (Z >= 0 && Z < NDENSELEM) ? indexZ[Z] : -1; } G4int G4DensityEffectData::GetIndex(const G4String& matName) const diff --git a/source/materials/src/G4IonisParamMat.cc b/source/materials/src/G4IonisParamMat.cc index 90fdce13f7..34d1795cda 100644 --- a/source/materials/src/G4IonisParamMat.cc +++ b/source/materials/src/G4IonisParamMat.cc @@ -212,11 +212,12 @@ void G4IonisParamMat::ComputeDensityEffectParameters() G4double corr = 0.0; if(idx < 0 && 1 == nelm) { - idx = fDensityData->GetElementIndex(Z0, fMaterial->GetState()); + G4int z = (1 == Z0 && State == kStateLiquid) ? 0 : Z0; + idx = fDensityData->GetElementIndex(z); // Correction for base material or for non-nominal density // Except cases of very different density defined in user code - if(idx >= 0) { + if(idx >= 0 && 0 < z) { G4double dens = nist->GetNominalDensity(Z0); if(dens <= 0.0) { idx = -1; } else { diff --git a/source/materials/src/G4Material.cc b/source/materials/src/G4Material.cc index d834476d96..7e3006322c 100644 --- a/source/materials/src/G4Material.cc +++ b/source/materials/src/G4Material.cc @@ -167,15 +167,6 @@ G4Material::G4Material(const G4String& name, G4double density, fPressure = pressure; fNbComponents = nComponents; - theElementVector = new G4ElementVector(); - theElementVector->reserve(fNbComponents); - - fAtomsVector = new G4int[fNbComponents]; - fMassFractionVector = new G4double[fNbComponents]; - for(G4int i=0; i; + fElm = new std::vector; + } if(fIdxComponent >= fNbComponents) { G4ExceptionDescription ed; ed << "For material " << fName << " and added element " @@ -381,51 +376,51 @@ G4Material::AddElementByNumberOfAtoms(const G4Element* elm, G4int nAtoms) } if(fMassFraction) { G4ExceptionDescription ed; - G4cout << "For material " << fName << " and added element " - << elm->GetName() << " with Natoms=" << nAtoms - << " problem: cannot add by number of atoms after " - << "addition of elements by mass fraction"; + ed << "For material " << fName << " and added element " + << elm->GetName() << " with Natoms=" << nAtoms + << " problem: cannot add by number of atoms after " + << "addition of elements by mass fraction"; G4Exception ("G4Material::AddElementByNumberOfAtoms()", "mat031", FatalException, ed, ""); } - G4Element* element = const_cast(elm); - // filling - if (fIdxComponent < fNbComponents) { - G4bool isAdded = false; + // filling + G4bool isAdded = false; + if(!fElm->empty()) { for (G4int i=0; iGetName() << ", Natoms=" << nAtoms - << ", fIdxComponent=" << fIdxComponent - << " problem: attempt to add the same element, which already is at idx=" - << i << " with the Natoms=" << fAtomsVector[i]; - G4Exception ("G4Material::AddElementByNumberOfAtoms()", "mat031", - JustWarning, ed, ""); - fAtomsVector[i] += nAtoms; + if ( elm == (*fElm)[i] ) { + (*fAtoms)[i] += nAtoms; + isAdded = true; break; } } - if(!isAdded) { - theElementVector->push_back(element); - fAtomsVector[fNumberOfElements] = nAtoms; - ++fNumberOfElements; - } + } + if(!isAdded) { + fElm->push_back(elm); + fAtoms->push_back(nAtoms); + ++fNumberOfElements; } ++fIdxComponent; - // is filled + + // is filled - complete composition of atoms if (fIdxComponent == fNbComponents) { - // compute proportion by mass + theElementVector = new G4ElementVector(); + theElementVector->reserve(fNumberOfElements); + fAtomsVector = new G4int[fNumberOfElements]; + fMassFractionVector = new G4double[fNumberOfElements]; + G4double Amol = 0.; for (G4int i=0; iGetA(); + theElementVector->push_back((*fElm)[i]); + fAtomsVector[i] = (*fAtoms)[i]; + G4double w = fAtomsVector[i]*(*fElm)[i]->GetA(); Amol += w; fMassFractionVector[i] = w; } for (G4int i=0; iGetName() << ", massFraction= " << fraction << ", fIdxComponent=" << fIdxComponent - << " problem: attempt to add more than the declared number of elements " + << "; attempt to add more than the declared number of components " << fIdxComponent << " >= " << fNbComponents; G4Exception ("G4Material::AddElementByMassFraction()", "mat031", FatalException, ed, ""); } - G4Element* element = const_cast(elm); + if(0 == fIdxComponent) { + fElmFrac = new std::vector; + fElm = new std::vector; + } // filling - if (fIdxComponent < fNbComponents) { - G4bool isAdded = false; + G4bool isAdded = false; + if(!fElm->empty()) { for (G4int i=0; iGetName() << ", massFraction= " << fraction - << ", fIdxComponent=" << fIdxComponent - << " problem: attempt to add the same element, which is already at idx=" - << i << " with the fraction " << fMassFractionVector[i]; - G4Exception ("G4Material::AddElementByMassFraction()", "mat031", - JustWarning, ed, ""); - fMassFractionVector[i] += fraction; + if ( elm == (*fElm)[i] ) { + (*fElmFrac)[i] += fraction; isAdded = true; break; } } - if(!isAdded) { - theElementVector->push_back(element); - fMassFractionVector[fNumberOfElements] = fraction; - ++fNumberOfElements; - } + } + if(!isAdded) { + fElm->push_back(elm); + fElmFrac->push_back(fraction); + ++fNumberOfElements; } ++fIdxComponent; @@ -527,15 +517,40 @@ void G4Material::AddMaterial(G4Material* material, G4double fraction) G4ExceptionDescription ed; ed << "For material " << fName << " and added material " << material->GetName() << ", massFraction= " << fraction - << " attempt to add more than the declared number of elements " + << "; attempt to add more than the declared number of components " << fIdxComponent << " >= " << fNbComponents; G4Exception ("G4Material::AddMaterial()", "mat031", FatalException, ed, ""); } + + if(0 == fIdxComponent) { + fElmFrac = new std::vector; + fElm = new std::vector; + } + // filling - if (fIdxComponent < fNbComponents) { - fMatComponents[material] = fraction; + G4int nelm = material->GetNumberOfElements(); + for(G4int j=0; jGetElement(j); + auto frac = material->GetFractionVector(); + G4bool isAdded = false; + if(!fElm->empty()) { + for (G4int i=0; ipush_back(elm); + fElmFrac->push_back(fraction*frac[j]); + ++fNumberOfElements; + } } + + fMatComponents[material] = fraction; ++fIdxComponent; // is filled @@ -547,70 +562,20 @@ void G4Material::AddMaterial(G4Material* material, G4double fraction) void G4Material::FillVectors() { // there are material components - if(!fMatComponents.empty()) { - G4int nel = fNumberOfElements; - // check list of materials - for(auto & x : fMatComponents) { - const G4Material* mat = x.first; - G4int nn = mat->GetNumberOfElements(); - for(G4int j=0; jGetElement(j); - for(G4int k=0; k fNbComponents) { - delete [] fAtomsVector; - fAtomsVector = new G4int[nel]; - G4double* v = new G4double[nel]; - for(G4int i=0; iGetNumberOfElements(); - const G4double* elmFrac = mat->GetFractionVector(); - for(G4int j=0; jGetElement(j); - for(G4int k=0; kpush_back(const_cast(elm)); - fMassFractionVector[fNumberOfElements] = frac*elmFrac[j]; - ++fNumberOfElements; - } - } - } - } - - // check sum of weights -- OK? + theElementVector = new G4ElementVector(); + theElementVector->reserve(fNumberOfElements); + fAtomsVector = new G4int[fNumberOfElements]; + fMassFractionVector = new G4double[fNumberOfElements]; G4double wtSum(0.0); for (G4int i=0; ipush_back((*fElm)[i]); + fMassFractionVector[i] = (*fElmFrac)[i]; wtSum += fMassFractionVector[i]; } + delete fElmFrac; + delete fElm; + + // check sum of weights -- OK? if (std::abs(1.-wtSum) > perThousand) { G4ExceptionDescription ed; ed << "For material " << fName << " sum of fractional masses " diff --git a/source/materials/src/G4NistMaterialBuilder.cc b/source/materials/src/G4NistMaterialBuilder.cc index d90177b2af..97c0d48d37 100644 --- a/source/materials/src/G4NistMaterialBuilder.cc +++ b/source/materials/src/G4NistMaterialBuilder.cc @@ -734,7 +734,7 @@ void G4NistMaterialBuilder::NistSimpleMaterials() AddMaterial("G4_Li", 0.534 , 3, 40. ); AddMaterial("G4_Be", 1.848 , 4, 63.7); AddMaterial("G4_B" , 2.37 , 5, 76. ); - AddMaterial("G4_C" , 2. , 6, 81. ); + AddMaterial("G4_C" , 2. , 6, 78. ); AddMaterial("G4_N" , 1.16520e-3, 7, 82. , 1, kStateGas); AddMaterial("G4_O" , 1.33151e-3, 8, 95. , 1, kStateGas); AddMaterial("G4_F" , 1.58029e-3, 9, 115. , 1, kStateGas);