voz1b1_2 compiles and links now.

This commit is contained in:
Guilhem Lavaux 2013-06-11 10:26:50 +02:00
parent 51efae2e09
commit 26132586c0
3 changed files with 72 additions and 52 deletions

View file

@ -26,6 +26,7 @@ bool checkParameters(int *numdiv, int *b)
return true;
}
struct BoxData
{
float width[3], totwidth[3];
@ -40,6 +41,8 @@ struct BoxData
void prepareBox(const PositionData& pdata, int *b_id, int *numdivs);
void checkParticle(double *xyz, bool& in_main, bool& in_buf);
void addGuardPoints();
void findBoundary();
};
void BoxData::checkParticle(double *xyz, bool& in_main, bool& in_buf)
@ -189,9 +192,9 @@ void BoxData::addGuardPoints()
{
/* Add guard points */
for (i=0; i<=NGUARD; i++)
for (int i=0; i <= NGUARD; i++)
{
for (j=0; j<=NGUARD; j++)
for (int j=0; j <= NGUARD; j++)
{
/* Bottom */
for (int a = 0; a < 3; a++)
@ -209,6 +212,52 @@ void BoxData::addGuardPoints()
}
}
void BoxData::findBoundary()
{
double BF = std::numeric_limits<double>::max();
for (int j = 0; j < 3; j++)
{
xyz_min[j] = -BF;
xyz_max[j] = BF;
}
for (pid_t i = nvpbuf; i < nvpall; i++) {
for (int j = 0; j < 3; j++)
{
xyz_min[j] = std::min(xyz_min[j], parts[3*i+j]);
xyz_max[j] = std::min(xyz_max[j], parts[3*i+j]);
}
}
}
void saveTesselation(const string& outfile, PositionData& pdata, BoxData& boxdata, PARTADJ *adjs, float *vols)
{
ofstream out(outfile.c_str());
if (!out)
{
cout << format("Unable to open %s") % outfile << endl;
exit(0);
}
out.write((char *)&pdata.np, sizeof(pid_t));
out.write((char *)&boxdata.nvp, sizeof(pid_t));
cout << format("nvp = %d") % boxdata.nvp << endl;
/* Tell us where the original particles were */
out.write((char *)boxdata.orig, sizeof(pid_t)*boxdata.nvp);
/* Volumes*/
out.write((char *)vols,sizeof(float)*boxdata.nvp);
/* Adjacencies */
for (pid_t i = 0; i < boxdata.nvp; i++)
{
out.write((char*)&adjs[i].nadj, sizeof(pid_t));
if (adjs[i].nadj > 0)
out.write((char *)adjs[i].adj, adjs[i].nadj*sizeof(pid_t));
else
(cout << "0").flush();
}
out.close();
}
int main(int argc, char *argv[]) {
PositionData pdata;
@ -270,26 +319,18 @@ int main(int argc, char *argv[]) {
else
boxdata.bf = 0.1;
boxdata.prepareBox(pdata);
boxdata.prepareBox(pdata, b, numdiv);
pdata.destroy();
boxdata.addGuardPoints();
adjs = new PARTADJ[np];
adjs = new PARTADJ[boxdata.nvpall];
if (adjs == 0)
{
cout << "Unable to allocate adjs" << endl;
return 0;
}
xmin = BF; xmax = -BF; ymin = BF; ymax = -BF; zmin = BF; zmax = -BF;
for (i=nvpbuf;i<nvpall;i++) {
if (parts[3*i] < xmin) xmin = parts[3*i];
if (parts[3*i] > xmax) xmax = parts[3*i];
if (parts[3*i+1] < ymin) ymin = parts[3*i+1];
if (parts[3*i+1] > ymax) ymax = parts[3*i+1];
if (parts[3*i+2] < zmin) zmin = parts[3*i+2];
if (parts[3*i+2] > zmax) zmax = parts[3*i+2];
}
boxdata.findBoundary();
cout << format("Added guard points to total %d points (should be %d)")
% boxdata.nvpall % (boxdata.nvpbuf + 6*(NGUARD+1)*(NGUARD+1)) << endl;
@ -297,7 +338,7 @@ int main(int argc, char *argv[]) {
/* Do tesselation*/
printf("File read. Tessellating ...\n"); fflush(stdout);
exitcode = delaunadj(boxdata.parts, boxdata.nvp, boxdata.nvpbuf, boxdata.nvpall, &adjs);
int exitcode = delaunadj(boxdata.parts, boxdata.nvp, boxdata.nvpbuf, boxdata.nvpall, &adjs);
if (exitcode != 0)
{
printf("Error while tesselating. Stopping here."); fflush(stdout);
@ -306,9 +347,9 @@ int main(int argc, char *argv[]) {
/* Now calculate volumes*/
printf("Now finding volumes ...\n"); fflush(stdout);
vols = new float[nvp];
vols = new float[boxdata.nvp];
for (pid_t i = 0; i < nvp; i++)
for (pid_t i = 0; i < boxdata.nvp; i++)
{ /* Just the original particles
* Assign adjacency coordinate array*/
/* Volumes */
@ -316,7 +357,7 @@ int main(int argc, char *argv[]) {
{
for (int d = 0; d < 3; d++)
{
deladjs[3*j + d] = parts[3*adjs[i].adj[j]+d] - parts[3*i+d];
deladjs[3*j + d] = boxdata.parts[3*adjs[i].adj[j]+d] - boxdata.parts[3*i+d];
if (deladjs[3*j+d] < -boxdata.width2[d])
deladjs[3*j+d] += boxdata.width[d];
@ -326,51 +367,28 @@ int main(int argc, char *argv[]) {
}
exitcode = vorvol(deladjs, points, intpoints, adjs[i].nadj, &(vols[i]));
vols[i] *= np/V0;
vols[i] *= pdata.np/pdata.V0;
if ((i % 1000) == 0)
cout << format("%d: %d, %f") % i % adjs[i].nadj % vols[i] << endl;
}
/* Get the adjacencies back to their original values */
for (pid_t i=0; i<nvp; i++)
for (pid_t i=0; i<boxdata.nvp; i++)
for (int j=0; j < adjs[i].nadj; j++)
adjs[i].adj[j] = orig[adjs[i].adj[j]];
totalvol = 0.;
for (pid_t i=0;i<nvp; i++)
for (pid_t i=0;i<boxdata.nvp; i++)
totalvol += vols[i];
cout << format("Average volume = %g") % (totalvol/nvp) << endl;
cout << format("Average volume = %g") % (totalvol/boxdata.nvp) << endl;
/* Now the output!
First number of particles */
outfile = str(format("%s/part.%s.%02d.%02d.%02d") % outDir % suffix % b[0] % b[1] % b[2]);
cout << format("Output to %s") %outfile << endl << endl;
out.open(outfile.c_str());
if (!out)
{
cout << format("Unable to open %s") % outfile << endl;
return 0;
}
out.write((char *)&pdata.np, sizeof(pid_t));
out.write((char *)&boxdata.nvp, sizeof(pid_t));
cout << format("nvp = %d") % nvp << endl;
/* Tell us where the original particles were */
out.write((char *)boxdata.orig, sizeof(pid_t)*boxdata.nvp);
/* Volumes*/
out.write((char *)vols,sizeof(float)*nvp);
/* Adjacencies */
for (pid_t i = 0; i < nvp; i++)
{
out.write((char*)&adjs[i].nadj, sizeof(pid_t));
if (adjs[i].nadj > 0)
out.write((char *)adjs[i].adj, adjs[i].nadj*sizeof(pid_t));
else
(cout << "0").flush();
}
out.close();
saveTesselation(outfile, pdata, boxdata, adjs, vols);
delete[] adjs;
return(0);

View file

@ -1,3 +1,4 @@
#include <limits>
#include <iostream>
#include <fstream>
#include "voz_io.hpp"
@ -12,15 +13,15 @@ bool PositionData::readFrom(const string& fname)
try
{
UnformattedRead f(fname.c_str());
f.beginCheckpoint();
np = f.readInt32();
f.endCheckPoint();
f.endCheckpoint();
for (int j = 0; j < 3; j++)
{
xyz[j] = new float[np];
f.beginCheckpoint();
for (int p = 0; p < np; p++)
xyz[j][p] = f.readReal32();
@ -39,7 +40,7 @@ bool PositionData::readFrom(const string& fname)
{
return false;
}
return true;
}
@ -49,13 +50,13 @@ void PositionData::findExtrema()
for (int i = 0; i < 3; i++)
{
xyz_min[i] = BF;
xyz_min[i] = BF;
xyz_max[i] = -BF;
}
for (int i = 0; i < 3; i++)
{
for (pid_t p = 0; p < np; p++)
for (pid_t p = 0; p < np; p++)
{
xyz_min[p] = min(xyz_min[p], xyz[p][i]);
xyz_max[p] = max(xyz_max[p], xyz[p][i]);

View file

@ -8,7 +8,8 @@ struct PositionData
float *xyz[3];
pid_t np;
float xyz_min[3], xyz_max[3];
float V0;
void destroy()
{
for (int j = 0; j < 3; j++)