Commit d2e81d3c authored by gsell's avatar gsell

Merge branch 'master' of gitorious.psi.ch:opal/src

parents 8792b0ca 03fb919c
......@@ -126,4 +126,8 @@ src/ValueDefinitions/cmake_install.cmake
src/cmake_install.cmake
src/config.h
src/opal
*~
\ No newline at end of file
unit_tests/gtest/include/
unit_tests/gtest/lib/
unit_tests/gtest/src/
*~
/Debug/
\ No newline at end of file
......@@ -19,8 +19,8 @@
//
// ------------------------------------------------------------------------
#include "Algebra/Tps.cpp"
#include "Algebra/Vps.cpp"
#include <Algebra/Tps.hpp>
#include <Algebra/Vps.hpp>
#include <complex>
......
......@@ -19,8 +19,8 @@
//
// ------------------------------------------------------------------------
#include "Algebra/Tps.cpp"
#include "Algebra/Vps.cpp"
#include <Algebra/Tps.hpp>
#include <Algebra/Vps.hpp>
// Force instantiation of Tps<double> class.
......
......@@ -19,8 +19,8 @@
//
// ------------------------------------------------------------------------
#include "Algebra/Tps.cpp"
#include "Algebra/Vps.cpp"
#include <Algebra/Tps.hpp>
#include <Algebra/Vps.hpp>
// Force instantiation of Tps< Tps<T> > class.
......
......@@ -209,6 +209,6 @@ PoissonBracket(const FLieGenerator<T, N> &, const FLieGenerator<T, N> &);
template <class T, int N>
std::ostream &operator<<(std::ostream &, const FLieGenerator<T, N> &);
#include "FixedAlgebra/FLieGenerator.cpp"
#include <FixedAlgebra/FLieGenerator.hpp>
#endif // CLASSIC_FLieGenerator_HH
......@@ -486,6 +486,6 @@ std::ostream &operator<<(std::ostream &os, const FTps<T, N> &);
// Implementation.
#include "FixedAlgebra/FVps.h"
#include "FixedAlgebra/FTps.cpp"
#include <FixedAlgebra/FTps.hpp>
#endif // CLASSIC_FTps_HH
......@@ -318,6 +318,6 @@ std::ostream &operator<<(std::ostream &os, const FVps<T, N> &vps);
// Implementation.
#include "FixedAlgebra/FVps.cpp"
#include "FixedAlgebra/FVps.hpp"
#endif // CLASSIC_FVps_HH
......@@ -733,9 +733,10 @@ FVps<T, N> FVps<T, N>::substitute(const FMatrix<T, N, N> &mat, int nl, int nh) c
// Initialisations.
// Array element fp[k][m] points to the next order m monomial
// to transform in k-th component.
const T *fp[N][nh+1];
const T **fp[N];
T *g[N];
for(int k = N; k-- > 0;) {
fp[k] = new const T* [nh+1];
for(int m = nl; m <= nh; ++m) fp[k][m] = f[k].begin(m);
g[k] = result[k].begin();
}
......@@ -831,6 +832,10 @@ FVps<T, N> FVps<T, N>::substitute(const FMatrix<T, N, N> &mat, int nl, int nh) c
oldvrbl = vrbl;
}
for(int k = N; k-- > 0;) {
delete[] fp[k];
}
return result;
}
......@@ -875,14 +880,16 @@ FVps<T, N> FVps<T, N>::substitute(const FVps<T, N> &rhs, int trunc) const {
if(nl == 0) nl = 1;
// Allocate working arrays.
const T *fp[N][nh+1];
const T ** fp[N];
Array1D< FTps<T, N> > t(nh + 1);
// Initialisations.
// Array element fp[k][m] points to the next order m monomial
// to transform in the k-th component.
for(int k = N; k-- > 0;)
for(int k = N; k-- > 0;) {
fp[k] = new const T* [nh+1];
for(int m = nl; m <= nh; ++m) fp[k][m] = f[k].begin(m);
}
const Array1D<int> *oldvrbl = 0;
int start_nh = FTps<T, N>::orderStart(nh), end_nh = FTps<T, N>::orderEnd(nh);
int nh1 = nh - 1, nh2 = nh - 2;
......@@ -963,6 +970,10 @@ FVps<T, N> FVps<T, N>::substitute(const FVps<T, N> &rhs, int trunc) const {
oldvrbl = vrbl;
}
for(int k = N; k-- > 0;) {
delete[] fp[k];
}
return result;
}
......
......@@ -234,7 +234,7 @@ std::ostream &operator<<(std::ostream &os, const LinearFun<T, N> &);
// Implementation.
#include "FixedAlgebra/LinearFun.cpp"
#include <FixedAlgebra/LinearFun.hpp>
#include "FixedAlgebra/LinearMap.h"
#endif // CLASSIC_LinearFun_HH
......@@ -214,6 +214,6 @@ std::ostream &operator<<(std::ostream &os, const LinearMap<T, N> &vps);
// Implementation.
#include "FixedAlgebra/LinearMap.cpp"
#include <FixedAlgebra/LinearMap.hpp>
#endif // CLASSIC_LinearMap_HH
......@@ -147,6 +147,6 @@ PoissonBracket(const Taylor<T> &, const Taylor<T> &);
template <class T>
std::ostream &operator<<(std::ostream &, const Taylor<T> &);
#include "FixedAlgebra/Taylor.cpp"
#include <FixedAlgebra/Taylor.hpp>
#endif // CLASSIC_Taylor_HH
......@@ -262,7 +262,7 @@ std::ostream &operator<<(std::ostream &os, const TransportFun<T, N> &);
// Implementation.
#include "FixedAlgebra/TransportFun.cpp"
#include <FixedAlgebra/TransportFun.hpp>
#include "FixedAlgebra/TransportMap.h"
#endif // CLASSIC_TransportFun_HH
......@@ -217,6 +217,6 @@ operator-(const FVector<T, N> &lhs, const TransportMap<T, N> &rhs);
// Implementation.
#include "FixedAlgebra/TransportMap.cpp"
#include <FixedAlgebra/TransportMap.hpp>
#endif // CLASSIC_TransportMap_HH
......@@ -7,10 +7,13 @@
/*****************************************************************************/
#include <algorithm>
#include <vector>
#include <cstring>
#include "Ctunes.h"
#include "lomb.h"
#include "lomb.cc"
#include "Utility/Inform.h"
#include "Algorithms/Ctunes.h"
#include "Algorithms/lomb.h"
#include "Algorithms/lomb.hpp"
extern Inform *gmsg;
......@@ -18,7 +21,11 @@ extern Inform *gmsg;
//RANLIB_class rndm(265314159,4);
TUNE_class::TUNE_class()
TUNE_class::TUNE_class():
ofac(0.0),
hifac(0.0),
Qmin(0.0),
Qmax(0.0)
/*---------------------------------------------------------------------------*
* constructor
* ===========
......@@ -101,7 +108,7 @@ int TUNE_class::LombAnalysis(vector<double> &x, vector<double> &y, int nhis, dou
memset(mess, '\0', sizeof(mess));
sprintf(mess, "@C3ERROR: Lomb analysis failed!");
*gmsg << "* " << mess << endl;
delete la;
la = NULL;
return(-1);
......@@ -240,7 +247,7 @@ int TUNE_class::LombAnalysis(double *x, double *y, int Ndat, int nhis)
memset(mess, '\0', sizeof(mess));
sprintf(mess, "%12.8f %8.2f %8.3f %d", pairx[pairc], pairy[pairc],
probi, i);
*gmsg << "* " << mess << endl;
*gmsg << "* " << mess << endl;
}
}
pairc++;
......@@ -256,5 +263,4 @@ int TUNE_class::LombAnalysis(double *x, double *y, int Ndat, int nhis)
return(0);
}
}
\ No newline at end of file
......@@ -16,6 +16,7 @@
//
// ------------------------------------------------------------------------
#include <Algorithms/Ctunes.hpp>
#include "Algorithms/ParallelCyclotronTracker.h"
#include <cfloat>
......@@ -71,7 +72,6 @@
#include "Utilities/OpalOptions.h"
#include "Ctunes.h"
#include "Ctunes.cc"
#include <cassert>
#include <hdf5.h>
......
This diff is collapsed.
......@@ -124,13 +124,17 @@ public:
// If [b]revBeam[/b] is true, the beam runs from s = C to s = 0.
// If [b]revTrack[/b] is true, we track against the beam.
explicit ParallelTTracker(const Beamline &bl, PartBunch &bunch, DataSink &ds,
const PartData &data, bool revBeam, bool revTrack, int maxSTEPS, double zstop, int timeIntegrator, size_t N);
const PartData &data, bool revBeam, bool revTrack,
const std::vector<unsigned long long> &maxSTEPS, const std::vector<double> &zstop,
int timeIntegrator, const std::vector<double> &dt, size_t N);
/// Constructor
// Amr pointer is taken
#ifdef HAVE_AMR_SOLVER
explicit ParallelTTracker(const Beamline &bl, PartBunch &bunch, DataSink &ds,
const PartData &data, bool revBeam, bool revTrack, int maxSTEPS, double zstop, int timeIntegrator, size_t N, Amr* amrptr_in);
const PartData &data, bool revBeam, bool revTrack,
const std::vector<unsigned long long> &maxSTEPS, const std::vector<double> &zstop,
int timeIntegrator, const std::vector<double> &dt, size_t N, Amr* amrptr_in);
#endif
virtual ~ParallelTTracker();
......@@ -282,7 +286,7 @@ private:
bool nEmissionMode_m;
/// where to stop
double zStop_m;
std::vector<double> zStop_m;
/// The scale factor for dimensionless variables (FIXME: move to PartBunch)
double scaleFactor_m;
......@@ -294,7 +298,8 @@ private:
double rescale_coeff_m;
double dtTrack_m;
double dtCurrentTrack_m;
std::vector<double> dtAllTracks_m;
double surfaceEmissionStop_m;
......@@ -324,7 +329,7 @@ private:
unsigned int emissionSteps_m;
/// The maximal number of steps the system is integrated per TRACK
unsigned long long localTrackSteps_m;
std::vector<unsigned long long> localTrackSteps_m;
size_t maxNparts_m;
size_t numberOfFieldEmittedParticles_m;
......@@ -419,11 +424,11 @@ private:
void setTime();
void initializeBoundaryGeometry();
void doBinaryRepartition();
void Tracker_Default();
void executeDefaultTracker();
#ifdef HAVE_AMR_SOLVER
void Tracker_AMR();
void executeAMRTracker();
#endif
void Tracker_AMTS();
void executeAMTSTracker();
void push(double h);
void kick(double h, bool avoidGammaCalc = false);
void computeExternalFields_AMTS();
......@@ -835,8 +840,8 @@ inline void ParallelTTracker::writePhaseSpace(const long long step, const double
msg << "* Wrote beam statistics." << endl;
itsDataSink_m->writeStatData(*itsBunch, FDext, rmax(2), sposRef, rmin(2), collimatorLosses);
}
// itsBunch->printBinHist();
}
#endif // OPAL_ParallelTTracker_HH
#endif // OPAL_ParallelTTracker_HH
\ No newline at end of file
......@@ -8,6 +8,7 @@
/*****************************************************************************/
#include <vector>
#include <cmath>
#include <functional>
typedef struct {
double x, y;
......
......@@ -6,6 +6,7 @@
/* (update: ASM, September 2001) */
/*****************************************************************************/
#include "lomb.h"
#include <iostream>
using namespace std;
......
......@@ -276,7 +276,7 @@ int main(int argc, char *argv[]) {
*gmsg << endl << "*** User error detected by function \""
<< ex.where() << "\"" << endl;
abort();
} catch(std::bad_alloc) {
} catch(std::bad_alloc &) {
*gmsg << "Sorry, virtual memory exhausted." << endl;
abort();
} catch(...) {
......
......@@ -590,6 +590,11 @@ void OpalParser::run() const {
*gmsg << " ";
stat->print();
*gmsg << " Sorry, virtual memory exhausted.\n" << endl;
} catch(assertion &ex) {
ERRORMSG("\n*** Runtime-error ******************\n");
ERRORMSG(ex.what());
ERRORMSG("\n************************************\n" << endl);
abort();
} catch(exception &ex) {
*gmsg << "\n*** Error:\n";
stat->printWhere(false);
......@@ -620,4 +625,4 @@ void OpalParser::run(TokenStream *is) const {
void OpalParser::stop() const {
stopFlag = true;
}
}
\ No newline at end of file
......@@ -2260,8 +2260,8 @@ void DataSink::writeSDDSHeader(ofstream &outputFile,
}
for (size_t i = 0; i < losses.size(); ++ i) {
outputFile << "&column name=" << losses[i].first << " losses, type=long, units=1, ";
outputFile << "description=\"" << columnStart++ << " " << losses[i].second << " losses \" &end" << endl;
outputFile << "&column name=" << losses[i].first << ", type=long, units=1, ";
outputFile << "description=\"" << columnStart++ << " " << losses[i].second << "\" &end" << endl;
}
outputFile << "&data mode=ascii &end" << endl;
......@@ -2895,4 +2895,4 @@ void DataSink::setOPALt() {
/***************************************************************************
* $RCSfile: DataSink.cpp,v $ $Author: adelmann $
* $Revision: 1.3 $ $Date: 2004/06/02 19:38:54 $
***************************************************************************/
\ No newline at end of file
***************************************************************************/
......@@ -36,8 +36,10 @@ otherwise a new bunch is allocated in the dictionary.
*/
Track::Track(BeamSequence *u, const PartData &ref, double dt, int maxtsteps, int stepsperturn, double zStop, int timeintegrator,
int nslices, double t0, double dtScInit, double deltaTau):
Track::Track(BeamSequence *u, const PartData &ref, const std::vector<double> & dt,
const std::vector<unsigned long long> & maxtsteps, int stepsperturn,
const std::vector<double> & zStop, int timeintegrator, int nslices,
double t0, double dtScInit, double deltaTau):
reference(ref),
use(u),
parser(),
......@@ -69,4 +71,4 @@ Track::Track(BeamSequence *u, const PartData &ref, double dt, int maxtsteps, int
Track::~Track()
{}
{}
\ No newline at end of file
......@@ -35,7 +35,9 @@ class Track {
public:
Track(BeamSequence *, const PartData &, double dt, int maxtsteps, int stepsperturn, double zStop, int timeintegrator, int nslices,
Track(BeamSequence *, const PartData &, const std::vector<double> & dt,
const std::vector<unsigned long long> & maxtsteps, int stepsperturn,
const std::vector<double> & zStop, int timeintegrator, int nslices,
double t0, double dtScInit, double deltaTau);
~Track();
......@@ -57,7 +59,7 @@ public:
static Track *block;
/// The initial timestep
double dT;
std::vector<double> dT;
// For AMTS integrator in OPAL-T
double dtScInit, deltaTau;
......@@ -67,13 +69,13 @@ public:
double t0_m;
/// Maximal number of timesteps
int localTimeSteps;
std::vector<unsigned long long> localTimeSteps;
/// The timsteps per revolution period. ONLY available for OPAL-cycl.
int stepsPerTurn;
/// The location at which the simulation stops
double zstop;
std::vector<double> zstop;
/// The ID of time integrator
// 0 --- RK-4(default)
......@@ -90,4 +92,4 @@ private:
void operator=(const Track &);
};
#endif // OPAL_Track_HH
#endif // OPAL_Track_HH
\ No newline at end of file
......@@ -59,20 +59,20 @@ TrackCmd::TrackCmd():
("LINE", "Name of lattice to be tracked");
itsAttr[BEAM] = Attributes::makeString
("BEAM", "Name of beam to be used", "UNNAMED_BEAM");
itsAttr[DT] = Attributes::makeReal
("DT", "THE INTEGRATION TIMESTEP IN SECONDS", 1e-12);
itsAttr[DT] = Attributes::makeRealArray
("DT", "THE INTEGRATION TIMESTEP IN SECONDS");
itsAttr[DTSCINIT] = Attributes::makeReal
("DTSCINIT", "Only for adaptive integrator: Initial time step for space charge integration", 1e-12);
itsAttr[DTAU] = Attributes::makeReal
("DTAU", "Only for adaptive integrator: Alternative way to set accuracy of space integration.", -1.0);
itsAttr[T0] = Attributes::makeReal
("T0", "THE ELAPSED TIME OF THE BUNCH IN SECONDS", 0.0);
itsAttr[MAXSTEPS] = Attributes::makeReal
("MAXSTEPS", "THE MAXIMUM NUMBER OF INTEGRATION STEPS DT, should be larger ZSTOP/(beta*c average)", 10);
itsAttr[MAXSTEPS] = Attributes::makeRealArray
("MAXSTEPS", "THE MAXIMUM NUMBER OF INTEGRATION STEPS DT, should be larger ZSTOP/(beta*c average)");
itsAttr[STEPSPERTURN] = Attributes::makeReal
("STEPSPERTURN", "THE TIME STEPS PER REVOLUTION PERIOD, ONLY FOR OPAL-CYCL", 720);
itsAttr[ZSTOP] = Attributes::makeReal
("ZSTOP", "Defines a z-location [m], after which the simulation stops when the last particles passes", 1000000.0);
itsAttr[ZSTOP] = Attributes::makeRealArray
("ZSTOP", "Defines a z-location [m], after which the simulation stops when the last particles passes");
itsAttr[TIMEINTEGRATOR] = Attributes::makeString
("TIMEINTEGRATOR", "Name of time integrator to be used", "RK-4");
itsAttr[NNB] = Attributes::makeReal
......@@ -92,8 +92,12 @@ TrackCmd *TrackCmd::clone(const std::string &name) {
return new TrackCmd(name, this);
}
double TrackCmd::getDT() const {
return Attributes::getReal(itsAttr[DT]);
std::vector<double> TrackCmd::getDT() const {
std::vector<double> dt = Attributes::getRealArray(itsAttr[DT]);
if (dt.size() == 0) {
dt.push_back(1e-12);
}
return dt;
}
double TrackCmd::getDTSCINIT() const {
......@@ -108,12 +112,30 @@ double TrackCmd::getT0() const {
return Attributes::getReal(itsAttr[T0]);
}
double TrackCmd::getZSTOP() const {
return Attributes::getReal(itsAttr[ZSTOP]);
std::vector<double> TrackCmd::getZSTOP() const {
std::vector<double> zstop = Attributes::getRealArray(itsAttr[ZSTOP]);
if (zstop.size() == 0) {
zstop.push_back(1000000.0);
}
return zstop;
}
int TrackCmd::getMAXSTEPS() const {
return (int) Attributes::getReal(itsAttr[MAXSTEPS]);
std::vector<unsigned long long> TrackCmd::getMAXSTEPS() const {
std::vector<double> maxsteps_d = Attributes::getRealArray(itsAttr[MAXSTEPS]);
std::vector<unsigned long long> maxsteps_i;
if (maxsteps_d.size() == 0) {
maxsteps_i.push_back(10ul);
}
for (auto it = maxsteps_d.begin(); it != maxsteps_d.end(); ++ it) {
if (*it < 0) {
maxsteps_i.push_back(10);
} else {
unsigned long long value = *it;
maxsteps_i.push_back(value);
}
}
return maxsteps_i;
}
int TrackCmd::getSTEPSPERTURN() const {
......@@ -147,19 +169,32 @@ void TrackCmd::execute() {
BeamSequence *use = BeamSequence::find(Attributes::getString(itsAttr[LINE]));
Beam *beam = Beam::find(Attributes::getString(itsAttr[BEAM]));
double dt = getDT();
std::vector<double> dt = getDT();
double t0 = getT0();
int maxsteps = getMAXSTEPS();
std::vector<unsigned long long> maxsteps = getMAXSTEPS();
int stepsperturn = getSTEPSPERTURN();
double zstop = getZSTOP();
std::vector<double> zstop = getZSTOP();
int timeintegrator = getTIMEINTEGRATOR();
int nslices = beam->getNumberOfSlices();
// Execute track block.
size_t numTracks = dt.size();
numTracks = std::max(numTracks, maxsteps.size());
numTracks = std::max(numTracks, zstop.size());
for (size_t i = dt.size(); i < numTracks; ++ i) {
dt.push_back(dt.back());
}
for (size_t i = maxsteps.size(); i < numTracks; ++ i) {
maxsteps.push_back(maxsteps.back());
}
for (size_t i = zstop.size(); i < numTracks; ++ i) {
zstop.push_back(zstop.back());
}
// Execute track block.
Track::block = new Track(use, beam->getReference(), dt, maxsteps, stepsperturn, zstop, timeintegrator, nslices, t0, getDTSCINIT(), getDTAU());
Track::block->parser.run();
// Clean up.
delete Track::block;
Track::block = 0;
}
}
\ No newline at end of file
......@@ -41,7 +41,7 @@ public:
virtual void execute();
/// Return the timestep in seconds
double getDT() const;
std::vector<double> getDT() const;
double getDTSCINIT() const;
......@@ -51,14 +51,14 @@ public:
double getT0() const;
/// Return the maximum timsteps we integrate the system
int getMAXSTEPS() const;
std::vector<unsigned long long> getMAXSTEPS() const;
/// Return the timsteps per revolution period. ONLY available for OPAL-cycl.
/// In OPAL-cycl, timestep is calculated by STEPSPERTURN, rather than given in TRACK command.
int getSTEPSPERTURN() const;
/// location at which the simulation stops
double getZSTOP() const;
std::vector<double> getZSTOP() const;
/// return the name of time integrator
int getTIMEINTEGRATOR() const;
......@@ -77,4 +77,4 @@ private:
TrackCmd(const std::string &name, TrackCmd *parent);
};
#endif // OPAL_TrackCmd_HH
#endif // OPAL_TrackCmd_HH
\ No newline at end of file
This diff is collapsed.
......@@ -69,7 +69,13 @@ private:
// Clone constructor.
TrackRun(const std::string &name, TrackRun *parent);
double SetDistributionParallelT(Beam *beam);
void setupSliceTracker();
void setupTTracker();
void setupCyclotronTracker();
void setupFieldsolver();
double setDistributionParallelT(Beam *beam);
ParallelTTracker *setupForAutophase();
// Pointer to tracking algorithm.
......@@ -87,9 +93,10 @@ private:
static const std::string defaultDistribution;
#ifdef HAVE_AMR_SOLVER
Amr* setupAMRSolver();
std::vector<std::string> filterString(std::string str);
std::pair<Box,unsigned int> getBlGrids(std::string str);
#endif
};
#endif // OPAL_TrackRun_HH
\ No newline at end of file
#endif // OPAL_TrackRun_HH
......@@ -40,7 +40,7 @@ RealConstant::RealConstant():
OpalData *OPAL = OpalData::getInstance();
OPAL->create(new RealConstant("PI", this, Physics::pi));
OPAL->create(new RealConstant("TWOPI", this, Physics::two_pi));
OPAL->create(new RealConstant("RADDEG", this, 180.0 / Physics::pi / 180.0));
OPAL->create(new RealConstant("RADDEG", this, 180.0 / Physics::pi));
OPAL->create(new RealConstant("DEGRAD", this, Physics::pi / 180.0));
OPAL->create(new RealConstant("E", this, Physics::e));
......
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment