Skip to content

Commit

Permalink
Fixed bug
Browse files Browse the repository at this point in the history
  • Loading branch information
@mseeger
mseeger committed Mar 22, 2014
1 parent 58a4add commit 6481743
Show file tree
Hide file tree
Showing 4 changed files with 30 additions and 17 deletions.
2 changes: 1 addition & 1 deletion src/eptools/potentials/ContainerPotManager.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -50,7 +50,7 @@
for (i=off=0; i<num; i++) {
startPos[i]=off; sz=parr[i]->size();
off+=sz;
nprec=pmArr[i]->numArgumentGroup(EPScalarPotential::atypeBivarPrec);
nprec=parr[i]->numArgumentGroup(EPScalarPotential::atypeBivarPrec);
if (havePrec && nprec<sz)
throw InvalidParameterException(EXCEPT_MSG("'atypeBivarPrec' potentials must form suffix"));
havePrec=(nprec>0);
Expand Down
8 changes: 7 additions & 1 deletion src/eptools/potentials/PotManagerFactory.cc
View file
Original file line number Diff line number Diff line change
Expand Up @@ -25,6 +25,7 @@
double* pvecP=parVec.p();
int* shrdP=parShrd.p();
bool hasBVPrec=false;
//char debStr[128]; // DEBUG

if (numPot.size()!=numk || annObj.size()!=numk)
throw InvalidParameterException(EXCEPT_MSG(""));
Expand All @@ -33,6 +34,7 @@
throw InvalidParameterException(EXCEPT_MSG(""));
if (numk>1)
parr.changeRep(numk);
//printMsgStdout("A");
for (k=0; k<numk; k++) {
// Create 'DefaultPotManager' object
//sprintf(debStr,"k=%d",k); printMsgStdout(debStr); // DEBUG!
Expand All @@ -41,13 +43,13 @@
// NOTE: We can pass 'pvecP' for construction parameters (if any),
// without having to prepare a parameter vector or even knowing the
// size. These parameters must form the prefix
//cout << "ManFactory: k=" << k << endl; // DEBUG!
Handle<EPScalarPotential> epPot;
try {
epPot.changeRep(EPPotentialFactory::createDefault(pid,pvecP,annObj[k]));
} catch (StandardException ex) {
throw InvalidParameterException(EXCEPT_MSG("Cannot create potential object"));
}
//sprintf(debStr,"B: epPot=%d",epPot.p()); printMsgStdout(debStr);
int numConstPars=epPot->numConstPars();
int npar=epPot->numPars(); // Could be 0
if (numConstPars>0) {
Expand All @@ -67,6 +69,8 @@
// Potential has no parameters
shrdMsk.changeRep(0); pvecMsk.changeRep(0);
}
// Ensure that if there are potentials of group 'atypeBivarPrec', they
// form the suffix
atype=epPot->getArgumentGroup();
if (hasBVPrec && atype!=EPScalarPotential::atypeBivarPrec)
throw InvalidParameterException(EXCEPT_MSG(""));
Expand All @@ -75,12 +79,14 @@
// they do not copy 'parShrd', 'parVec' content!
PotentialManager* pmanP=
new DefaultPotManager(epPot,npot,pvecMsk,shrdMsk,false);
//printMsgStdout("C");
if (numk==1)
return pmanP; // Single 'DefaultPotManager'
else
parr[k].changeRep(pmanP);
}
MYASS(numk>1);
//printMsgStdout("D");

return new ContainerPotManager(parr);
}
Expand Down
5 changes: 5 additions & 0 deletions src/eptools/wrap/eptools_helper.cc
View file
Original file line number Diff line number Diff line change
Expand Up @@ -27,6 +27,11 @@ void createPotentialManager(W_IARRAY(potids),W_IARRAY(numpot),W_DARRAY(parvec),
W_MASKARRAY(parvec);
W_MASKARRAY(parshrd);
W_MASKARRAY(annobj);
// DEBUG:
//sprintf(W_ERRSTR,"potids=%d,numpot=%d,parvec=%d,parshrd=%d,annobj=%d",
// potidsA.size(),numpotA.size(),parvecA.size(),parshrdA.size(),
// annobjA.size());
//printMsgStdout(W_ERRSTR);
try {
potMan.changeRep(PotManagerFactory::create(potidsA,numpotA,parvecA,
parshrdA,annobjA));
Expand Down
32 changes: 17 additions & 15 deletions src/eptools/wrap/eptwrap_fact_sequpdates.cc
View file
Original file line number Diff line number Diff line change
Expand Up @@ -124,21 +124,23 @@ void eptwrap_fact_sequpdates(int ain,int aout,int n,int m,W_IARRAY(updjind),
Interval<int> ivM(0,m-1,IntVal::ivClosed,IntVal::ivClosed);
if (ivM.check(updjind,nupdjind)!=0)
W_RETERROR(1,"UPDJIND: Entries of out range");
printMsgStdout("Point 1");
//printMsgStdout("Point 1");
/* Potential manager */
Handle<PotentialManager> potMan;
createPotentialManager(W_ARR(pm_potids),W_ARR(pm_numpot),W_ARR(pm_parvec),
W_ARR(pm_parshrd),W_ARR(pm_annobj),potMan,W_ERRARGS);
W_ARR(pm_parshrd),W_ARR(pm_annobj),potMan,
W_ERRARGS);
//sprintf(W_ERRSTR,"potMan=%d",potMan.p()); printMsgStdout(W_ERRSTR);
if (potMan->size()!=m)
W_RETERROR(1,"PM_*: Potential manager has wrong size");
/* Representation of B */
printMsgStdout("Point 2");
//printMsgStdout("Point 2");
Handle<FactorizedEPRepresentation> epRepr;
createFactEPRepres(n,m,W_ARR(rp_rowind),W_ARR(rp_colind),W_ARR(rp_bvals),
W_ARR(rp_pi),W_ARR(rp_beta),epRepr,W_ERRARGS);
/* Variable marginals */
ArrayHandle<double> margpiA,margbetaA;
printMsgStdout("Point 3");
//printMsgStdout("Point 3");
W_CHKSIZE(margpi,n,"MARGPI");
W_CHKSIZE(margpi,n,"MARGBETA");
W_MASKARRAY(margpi);
Expand All @@ -153,7 +155,7 @@ void eptwrap_fact_sequpdates(int ain,int aout,int n,int m,W_IARRAY(updjind),
W_RETERROR(1,"DAMPFACT: Out of range");
if (ain>17) {
// Selective damping
printMsgStdout("Point 4");
//printMsgStdout("Point 4");
if (ain<20)
W_RETERROR(1,"Need all SD_XXX or none");
W_CHKSIZE(sd_numvalid,n,"SD_NUMVALID");
Expand All @@ -164,7 +166,7 @@ void eptwrap_fact_sequpdates(int ain,int aout,int n,int m,W_IARRAY(updjind),
W_MASKARRAY(sd_topind);
W_CHKSIZE(sd_topval,nsd_topind,"SD_TOPVAL");
W_MASKARRAY(sd_topval);
printMsgStdout("Point 5");
//printMsgStdout("Point 5");
if (ain>20) {
if (nsd_subind==0 || nsd_subind>m)
W_RETERROR(1,"SD_SUBIND: Wrong size");
Expand Down Expand Up @@ -203,12 +205,12 @@ void eptwrap_fact_sequpdates(int ain,int aout,int n,int m,W_IARRAY(updjind),
}
/* Create max_pi data structure (only if selective damping) */
Handle<FactEPMaximumPiValues> epMaxPi;
printMsgStdout("Point 6");
//printMsgStdout("Point 6");
if (sd_k>0) {
try {
sprintf(W_ERRSTR,"MEX: n=%d,K=%d,numvalid=%d,topind=%d,topval=%d",n,
sd_k,sd_numvalidA.size(),sd_topindA.size(),sd_topvalA.size());
printMsgStdout(W_ERRSTR);
//sprintf(W_ERRSTR,"MEX: n=%d,K=%d,numvalid=%d,topind=%d,topval=%d",n,
// sd_k,sd_numvalidA.size(),sd_topindA.size(),sd_topvalA.size());
//printMsgStdout(W_ERRSTR);
epMaxPi.changeRep(new FactEPMaximumPiValues(epRepr,sd_k,
sd_numvalidA,sd_topindA,
sd_topvalA,sd_subindA,
Expand All @@ -221,7 +223,7 @@ void eptwrap_fact_sequpdates(int ain,int aout,int n,int m,W_IARRAY(updjind),
}
/* Create EP driver */
Handle<FactorizedEPDriver> epDriver;
printMsgStdout("Point 7");
//printMsgStdout("Point 7");
try {
epDriver.changeRep(new FactorizedEPDriver(potMan,epRepr,margbetaA,
margpiA,piminthres,epMaxPi));
Expand All @@ -234,8 +236,8 @@ void eptwrap_fact_sequpdates(int ain,int aout,int n,int m,W_IARRAY(updjind),
/* Main loop over updates */
for (int i=0; i<nupdjind; i++) {
int j=updjind[i];
sprintf(W_ERRSTR,"i=%d, j=%d",i,j);
printMsgStdout(W_ERRSTR);
//sprintf(W_ERRSTR,"i=%d, j=%d",i,j);
//printMsgStdout(W_ERRSTR);
int irstat=epDriver->sequentialUpdate(j,dampfact,(delta!=0)?(delta+i):0,
(sd_dampfact!=0)?(sd_dampfact+i):0);
if (rstat!=0) rstat[i]=irstat;
Expand All @@ -244,13 +246,13 @@ void eptwrap_fact_sequpdates(int ain,int aout,int n,int m,W_IARRAY(updjind),
if (irstat!=FactorizedEPDriver::updSuccess && sd_dampfact!=0)
sd_dampfact[i]=1.0;
}
printMsgStdout("Point 9");
//printMsgStdout("Point 9");
if (sd_nupd!=0) {
int inrec;
epMaxPi->getStats(*sd_nupd,inrec);
if (sd_nrec!=0) *sd_nrec=inrec;
}
printMsgStdout("Point 10");
//printMsgStdout("Point 10");
W_RETOK;
} catch (StandardException ex) {
W_RETERROR_ARGS(1,"Caught LHOTSE exception: %s", ex.msg());
Expand Down

0 comments on commit 6481743

Please sign in to comment.