cuIBM
A GPU-based Immersed Boundary Method code
updateRHS1.cu
Go to the documentation of this file.
1 
10 #include "updateRHS1.h"
11 
12 
17 namespace kernels
18 {
19 
29 __global__
30 void updateRHS1(real *rhs1, int numUV, int *tags)
31 {
32  int I = blockIdx.x*blockDim.x + threadIdx.x;
33 
34  if(I>=numUV)
35  return;
36 
37  rhs1[I] = rhs1[I]*(tags[I]==-1);
38 } // updateRHS1
39 
40 
53 __global__
54 void updateRHS1X(real *rhs1, int nx, int ny, real dt, real *dx, int *tags, real *coeffs, real *uv)
55 {
56  int I = blockIdx.x*blockDim.x + threadIdx.x;
57  int i = I % (nx-1);
58 
59  if( I < (nx-1)*ny )
60  {
61  rhs1[I] = (tags[I]==-1)*rhs1[I]
62  + ((tags[I]!=-1)*((1.0-coeffs[I])*uv[I])) * 0.5*(dx[i+1]+dx[i])/dt;
63  }
64 } // updateRHS1X
65 
66 
79 __global__
80 void updateRHS1Y(real *rhs1, int nx, int ny, real dt, real *dy, int *tags, real *coeffs, real *uv)
81 {
82  int numU = (nx-1)*ny;
83  int I = blockIdx.x*blockDim.x + threadIdx.x + numU;
84  int j = (I-numU) / nx;
85 
86  if( I < numU + nx*(ny-1) )
87  {
88  rhs1[I] = (tags[I]==-1)*rhs1[I]
89  + ((tags[I]!=-1)*((1.0-coeffs[I])*uv[I])) * 0.5*(dy[j+1]+dy[j])/dt;
90  }
91 } // updateRHS1Y
92 
93 
107 __global__
108 void updateRHS1X(real *rhs1, int nx, int ny, real dt, real *dx, int *tags, real *coeffs, real *coeffs2, real *uv)
109 {
110  int I = blockIdx.x*blockDim.x + threadIdx.x;
111  int i = I % (nx-1);
112 
113  if( I < (nx-1)*ny )
114  {
115  rhs1[I] = (tags[I]==-1)*rhs1[I]
116  + ((tags[I]!=-1)*((1.0-coeffs[I]-coeffs2[I])*uv[I])) * 0.5*(dx[i+1]+dx[i])/dt;
117  }
118 } // updateRHS1X
119 
120 
134 __global__
135 void updateRHS1Y(real *rhs1, int nx, int ny, real dt, real *dy, int *tags, real *coeffs, real *coeffs2, real *uv)
136 {
137  int numU = (nx-1)*ny;
138  int I = blockIdx.x*blockDim.x + threadIdx.x + numU;
139  int j = (I-numU) / nx;
140 
141  if( I < numU + nx*(ny-1) )
142  {
143  rhs1[I] = (tags[I]==-1)*rhs1[I]
144  + ((tags[I]!=-1)*((1.0-coeffs[I]-coeffs2[I])*uv[I])) * 0.5*(dy[j+1]+dy[j])/dt;
145  }
146 } // updateRHS1Y
147 
148 } // End of namespace kernels
Declaration of the kernels to update the right hand-side of the intermediate velocity flux solver...
__global__ void updateRHS1(real *rhs1, int numUV, int *tags)
Update the RHS vector of the velocity system at forcing nodes.
Definition: updateRHS1.cu:30
double real
Is a float or a double depending on the machine precision.
Definition: types.h:116
Contains all the custom-written CUDA kernels.
__global__ void updateRHS1Y(real *rhs1, int nx, int ny, real dt, real *dy, int *tags, real *coeffs, real *uv)
Perform 1D linear interpolation at forcing points for the y-velocity.
Definition: updateRHS1.cu:80
__global__ void updateRHS1X(real *rhs1, int nx, int ny, real dt, real *dx, int *tags, real *coeffs, real *uv)
Perform 1D linear interpolation at forcing points for the x-velocity.
Definition: updateRHS1.cu:54