Mentions légales du service

Skip to content
Snippets Groups Projects
Commit 0569d1fc authored by Olivier COULAUD's avatar Olivier COULAUD
Browse files

Modifications for CEA report

parent 0b6ed9f6
No related branches found
No related tags found
No related merge requests found
...@@ -30,7 +30,7 @@ ...@@ -30,7 +30,7 @@
//! <b> General arguments:</b> //! <b> General arguments:</b>
//! \param -help (-h) to see the parameters available in this driver //! \param -help (-h) to see the parameters available in this driver
//! \param -fin name: file name to convert (with extension .fma (ascii) or bfma (binary) //! \param -fin name: file name to convert (with extension .fma (ascii) or bfma (binary)
//! \param -fdlpoly name file coming from a DLpoly simulation //! \param -fdlpoly name file coming from a DLpoly simulation
//! //!
//! \param -fout name: generic name for files (without extension) and save data //! \param -fout name: generic name for files (without extension) and save data
//! with following format in name.fma or name.bfma if -bin is set" //! with following format in name.fma or name.bfma if -bin is set"
...@@ -46,11 +46,13 @@ ...@@ -46,11 +46,13 @@
int main(int argc, char ** argv){ int main(int argc, char ** argv){
FHelpDescribeAndExit(argc, argv, FHelpDescribeAndExit(argc, argv,
"Driver to change the format of the input file. " "Driver to change the format of the input file. "
"fdlpoly is not supported for now.", "fdlpoly is not supported for now.",
FParameterDefinitions::InputFile, FParameterDefinitions::OutputFile, FParameterDefinitions::InputFile, FParameterDefinitions::OutputFile,
FParameterDefinitions::OutputVisuFile, FParameterDefinitions::OutputBinFormat); FParameterDefinitions::OutputVisuFile,FParameterDefinitions::FormatVisuFile,
FParameterDefinitions::OutputBinFormat);
FSize NbPoints; FSize NbPoints;
FReal * particles = nullptr ; FReal * particles = nullptr ;
......
...@@ -34,10 +34,11 @@ int compareTwoArrays(const std::string &tag, const int &nbParticles, const clas ...@@ -34,10 +34,11 @@ int compareTwoArrays(const std::string &tag, const int &nbParticles, const clas
if(array1[idxPart].getPosition() != array2[idxPart].getPosition() ){ if(array1[idxPart].getPosition() != array2[idxPart].getPosition() ){
std::cerr <<"Wrong positions " <<std::endl std::cerr <<"Wrong positions " <<std::endl
<< " P1: " <<array1[idxPart].getPosition()<<std::endl << " P1: " <<array1[idxPart].getPosition()<<std::endl
<< " P2: " <<array2[idxPart].getPosition()<<std::endl; << " P2: " <<array2[idxPart].getPosition()<<std::endl
<< " error " <<array1[idxPart].getPosition()-array2[idxPart].getPosition()<<std::endl;
} }
potentialDiff.add(array1[idxPart].getPotential(),array2[idxPart].getPotential()); potentialDiff.add(array1[idxPart].getPotential(),array2[idxPart].getPotential());
fx.add(array1[idxPart].getForces()[0],array1[idxPart].getForces()[0]); fx.add(array1[idxPart].getForces()[0],array2[idxPart].getForces()[0]);
fy.add(array1[idxPart].getForces()[1],array2[idxPart].getForces()[1]); fy.add(array1[idxPart].getForces()[1],array2[idxPart].getForces()[1]);
fz.add(array1[idxPart].getForces()[2],array2[idxPart].getForces()[2]); fz.add(array1[idxPart].getForces()[2],array2[idxPart].getForces()[2]);
energy1 += array1[idxPart].getPhysicalValue() *array1[idxPart].getPotential() ; energy1 += array1[idxPart].getPhysicalValue() *array1[idxPart].getPotential() ;
......
...@@ -41,10 +41,10 @@ ...@@ -41,10 +41,10 @@
int main(int argc, char ** argv){ int main(int argc, char ** argv){
const FParameterNames LocalParameterEwald{ // const FParameterNames LocalParameterEwald{
{"-ewaldfile2"}, // {"-ewaldfile2"},
"-ewaldfile2 name2 if name2 contains the result done by the ewald method for 1/r kernel" // "-ewaldfile2 name2 if name2 contains the result done by the ewald method for 1/r kernel"
}; // };
FHelpDescribeAndExit(argc, argv, FHelpDescribeAndExit(argc, argv,
"Driver to change the format of the input file.", "Driver to change the format of the input file.",
FParameterDefinitions::InputFileOne, FParameterDefinitions::InputFileTwow, FParameterDefinitions::InputFileOne, FParameterDefinitions::InputFileTwow,
......
...@@ -50,7 +50,8 @@ void genDistusage() { ...@@ -50,7 +50,8 @@ void genDistusage() {
<< " -fout fileout.{bfma, fma) output file name "<< std::endl << " -fout fileout.{bfma, fma) output file name "<< std::endl
<< " -Ewald2FMM to add the Electric Moment in the potential and the force " << std::endl << " -Ewald2FMM to add the Electric Moment in the potential and the force " << std::endl
<< " -energy Value of the energy (Ewald)" << std::endl << " -energy Value of the energy (Ewald)" << std::endl
<< " -FMM2Ewald to remove the Electric Moment in the potential and the force " << std::endl ; << " -FMM2Ewald to remove the Electric Moment in the potential and the force " << std::endl
<< " Scaling parameters: -noscale, -dlpoly, -stamp" <<std::endl;
std::cout << " -verbose : print index x y z fx fy fy Q and V" << std::endl; std::cout << " -verbose : print index x y z fx fy fy Q and V" << std::endl;
exit(-1); exit(-1);
} }
...@@ -60,14 +61,9 @@ int main(int argc, char ** argv){ ...@@ -60,14 +61,9 @@ int main(int argc, char ** argv){
bool useDLPOLY=false,useFMA=false,useSTAMP=false; bool useDLPOLY=false,useFMA=false,useSTAMP=false;
// //
FReal scaleEnergy, scaleForce ; FReal scaleEnergy, scaleForce ;
const FReal r4pie0 = FReal(138935.4835);
const double q0 = 1.6021892e-19; const double q0 = 1.6021892e-19;
const double e0 = 8.854187e-12; const double e0 = 8.854187e-12;
const double ang = 1.0e-10; const double ang = 1.0e-10;
const double Unsur4pie0 =8.98755179e+9; // 1.0/(4*FMath::FPi*e0); //8.98755179e+9 ;
double D0 = 1.0;
std::cout << "Unsur4pie0: " << Unsur4pie0 << " 8.98755179e+9 Diff " << Unsur4pie0-8.98755179e+9 <<std::endl;
if( argc == 1 ){ if( argc == 1 ){
genDistusage() ; genDistusage() ;
...@@ -88,17 +84,27 @@ int main(int argc, char ** argv){ ...@@ -88,17 +84,27 @@ int main(int argc, char ** argv){
scaleForce = 1.0 ; scaleForce = 1.0 ;
// bool scale = true ; // bool scale = true ;
if(FParameters::existParameter(argc, argv, "-dlpoly")){ if(FParameters::existParameter(argc, argv, "-dlpoly")){
const FReal r4pie0 = FReal(138935.4835);
// scale = false ; // scale = false ;
scaleEnergy = r4pie0 / 418.4 ; // kcal mol^{-1} scaleEnergy = r4pie0 / 418.4 ; // kcal mol^{-1}
scaleForce = -r4pie0 ; // 10 J mol^{-1} A^{-1} scaleForce = -r4pie0 ; // 10 J mol^{-1} A^{-1}
D0 = 1.0;
useDLPOLY = true ; useDLPOLY = true ;
} }
else if(FParameters::existParameter(argc, argv, "-stamp")){ else if(FParameters::existParameter(argc, argv, "-stamp")){
// used in Stamp
const double Unsur4pie0 =8.98755179e+9; // 1.0/(4*FMath::FPi*e0); //8.98755179e+9 ;
scaleEnergy = q0*q0*Unsur4pie0 /ang; scaleEnergy = q0*q0*Unsur4pie0 /ang;
scaleForce = -scaleEnergy/ang; scaleForce = -scaleEnergy/ang;
D0 = q0/ang/ang; useSTAMP = true;
useSTAMP = true; }
else if(FParameters::existParameter(argc, argv, "-noscale")){
scaleEnergy = 1.0;
scaleForce = -1.0;
}
else {
std::cout << "Parameter -dlpoly, -stamp, -noscale is not set" <<std::endl;
std::exit(-1);
} }
std::cout <<"Scaling factor " <<std::endl std::cout <<"Scaling factor " <<std::endl
<< " Energy factor " << scaleEnergy<<std::endl << " Energy factor " << scaleEnergy<<std::endl
...@@ -184,14 +190,7 @@ int main(int argc, char ** argv){ ...@@ -184,14 +190,7 @@ int main(int argc, char ** argv){
} }
boxsize[0] = boxsize[1]= boxsize[2] = loaderFMA->getBoxWidth(); boxsize[0] = boxsize[1]= boxsize[2] = loaderFMA->getBoxWidth();
// scaleEnergy /= boxsize[0] ;
// scaleForce = -scaleEnergy/boxsize[0] ;
// boxsize[0] = boxsize[1]= boxsize[2] = 1.0 ;
numberofParticles = loaderFMA->getNumberOfParticles() ; numberofParticles = loaderFMA->getNumberOfParticles() ;
std::cout << " numberofParticles "<< std::endl;
particlesIn = new FmaRWParticle<8,8>[numberofParticles]; particlesIn = new FmaRWParticle<8,8>[numberofParticles];
memset(particlesIn, 0, sizeof( FmaRWParticle<8,8>) * numberofParticles) ; memset(particlesIn, 0, sizeof( FmaRWParticle<8,8>) * numberofParticles) ;
loaderFMA->fillParticle(particlesIn,numberofParticles); loaderFMA->fillParticle(particlesIn,numberofParticles);
...@@ -200,6 +199,7 @@ int main(int argc, char ** argv){ ...@@ -200,6 +199,7 @@ int main(int argc, char ** argv){
counter.tic(); counter.tic();
FPoint electricMoment(0.0,0.0,0.0) ; FPoint electricMoment(0.0,0.0,0.0) ;
FPoint xmin,xmax ;
// const --> then shared // const --> then shared
// //
...@@ -247,7 +247,7 @@ int main(int argc, char ** argv){ ...@@ -247,7 +247,7 @@ int main(int argc, char ** argv){
std::cout << " -Ewald2FMM or -FMM2Ewald should be set"<<std::endl; std::cout << " -Ewald2FMM or -FMM2Ewald should be set"<<std::endl;
exit(EXIT_FAILURE); exit(EXIT_FAILURE);
} }
std::cout << "coeffCorrection: "<<coeffCorrection<<" " <<coeffCorrection*D0<< std::endl; std::cout << "coeffCorrection: "<<coeffCorrection<< std::endl;
std::cout << "preScaleEnergy: "<<preScaleEnergy<<std::endl; std::cout << "preScaleEnergy: "<<preScaleEnergy<<std::endl;
std::cout << "postScaleEnergy: "<<postScaleEnergy<<std::endl; std::cout << "postScaleEnergy: "<<postScaleEnergy<<std::endl;
std::cout << "preScaleForce: "<<preScaleForce<<std::endl; std::cout << "preScaleForce: "<<preScaleForce<<std::endl;
...@@ -259,7 +259,7 @@ int main(int argc, char ** argv){ ...@@ -259,7 +259,7 @@ int main(int argc, char ** argv){
+ particlesIn[idx].getPosition().getZ()*electricMoment.getZ() ; + particlesIn[idx].getPosition().getZ()*electricMoment.getZ() ;
// Potential rescaling ( Not checked) // Potential rescaling ( Not checked)
double P = particlesIn[idx].getPotential(); double P = particlesIn[idx].getPotential();
particlesIn[idx].setPotential( P + tmp*coeffCorrection); // particlesIn[idx].setPotential( P + tmp*coeffCorrection);
// //
// Force rescaling // Force rescaling
double * forces ; double * forces ;
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment