Skip to content

Commit

Permalink
Bugfix for when 100% are outside the domain
Browse files Browse the repository at this point in the history
  • Loading branch information
CyprienBosserelle committed Feb 19, 2017
1 parent 771213f commit 68e65d0
Show file tree
Hide file tree
Showing 3 changed files with 28 additions and 18 deletions.
6 changes: 3 additions & 3 deletions Flow_kernel.cu
Original file line number Diff line number Diff line change
Expand Up @@ -1215,7 +1215,7 @@ __global__ void eulervstep(int nx, int ny, DECNUM dx, DECNUM dt, DECNUM g, DECNU
vvi[tx][ty] = 0.0f;
viscv[i] = 0.0f;
}
if (ix > 0 && iy>0 && iy<ny)
if (ix > 0)// && iy>0 && iy<ny)
{
vv[i] = vvi[tx][ty];
}//vdvdy[i]=tauby;
Expand Down Expand Up @@ -1448,12 +1448,12 @@ __global__ void uuvvzslatbnd(int nx, int ny, DECNUM * uu, DECNUM * vv, DECNUM *z
if (iy == ny - 1)
{
uu[i] = uub[tx][ty];
vv[i] = 0.0f;// vvb2[tx][ty];
vv[i] = vvb[tx][ty];
zs[i] = zsb[tx][ty];
}
if (iy==ny-2)
{
vv[i]=vvb[tx][ty];// THis is to follow XBeach definition although I don't really agree with it
//vv[i]=vvb[tx][ty];// THis is to follow XBeach definition although I don't really agree with it
// It should be that vv(i,ny-1)=vv(i,ny-2) end of story
}
if (ix == 0)
Expand Down
4 changes: 4 additions & 0 deletions XBeachGPU.h
Original file line number Diff line number Diff line change
Expand Up @@ -38,6 +38,10 @@ using DECNUM = float;
typedef std::complex<double> Complex;
typedef std::valarray<Complex> CArray;

// This below is to built a 3d array that may be too large for a simple malloc array Malloc c
// instead we can design a more complex structure of vectors


//class template to create vecotr ofg output locations fro the model
class TSnode{
public:
Expand Down
36 changes: 21 additions & 15 deletions makjonswap.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -258,7 +258,7 @@ void GenWGnLBW(XBGPUParam Param, int nf, int ndir,double * HRfreq,double * HRdir
double fmax,Sfmax; // Should be in Param
//int nspr = 0; // Should be in Param (need test for nspr ==1)
double * binedgeleft, * binedgeright; // size of ntheta
double * zeta, *Ampzeta; //water elevation ny*ntheta*tslen
//double * zeta, *Ampzeta; //water elevation ny*ntheta*tslen
double *eta, *Amp; //water elevation integrated over directions ny*tslen
double *stdzeta; //size of ntheta
double *E_tdir; // tslen
Expand Down Expand Up @@ -621,7 +621,11 @@ void GenWGnLBW(XBGPUParam Param, int nf, int ndir,double * HRfreq,double * HRdir
//= 0.0 + 1.0*I;
TwoDee<std::complex<double>> CompFn(ny, tslen);

zeta = (double *)malloc(ny*Param.ntheta*tslen*sizeof(double));
std::valarray<double> zeta(ny*Param.ntheta*tslen);



//zeta = (double *)malloc(ny*Param.ntheta*tslen*sizeof(double));
//Ampzeta = (double *)malloc(ny*Param.ntheta*tslen*sizeof(double));
eta = (double *)malloc(ny*tslen*sizeof(double));
Amp = (double *)malloc(ny*tslen*sizeof(double));
Expand Down Expand Up @@ -936,24 +940,26 @@ void GenWGnLBW(XBGPUParam Param, int nf, int ndir,double * HRfreq,double * HRdir
}

stdeta = sqrt(stdeta / (tslen - 1));
for (int itheta = 0; itheta < Param.ntheta; itheta++)
if (stdeta > 0.0)
{
temp = 0.0;

for (int n = 0; n < tslen; n++)
for (int itheta = 0; itheta < Param.ntheta; itheta++)
{
temp = temp + zeta[j + itheta*ny + n*ny*Param.ntheta] * zeta[j + itheta*ny + n*ny*Param.ntheta];
temp = 0.0;

}
stdzeta[itheta] = sqrt(temp / (tslen - 1));
for (int n = 0; n < tslen; n++)
{
temp = temp + zeta[j + itheta*ny + n*ny*Param.ntheta] * zeta[j + itheta*ny + n*ny*Param.ntheta];

for (int n = 0; n < tslen; n++)
{
zeta[j + itheta*ny + n*ny*Param.ntheta] = Amp[j + n*ny] * stdzeta[itheta] / stdeta;
}
}
stdzeta[itheta] = sqrt(temp / (tslen - 1));

}
for (int n = 0; n < tslen; n++)
{
zeta[j + itheta*ny + n*ny*Param.ntheta] = Amp[j + n*ny] * stdzeta[itheta] / stdeta;
}

}
}
//Calculate energy and interpolate to the output time step
// (Maybe dtin==dtbc or maybe not)
for (int itheta = 0; itheta < Param.ntheta; itheta++)
Expand Down Expand Up @@ -1382,7 +1388,7 @@ void GenWGnLBW(XBGPUParam Param, int nf, int ndir,double * HRfreq,double * HRdir
free(WDindex);

free(tin);
free(zeta);
//free(zeta); now a vector that will autodestruct??
//free(Ampzeta);//Reused zeta
free(eta);
free(Amp);
Expand Down

0 comments on commit 68e65d0

Please sign in to comment.