rifki sadikin 6 years ago
commit 8b1bf8d29e
  1. 56
      gpulib/PoissonSolver3DGPU.cu

@ -11,6 +11,17 @@ __device__ __constant__ float d_tempRatioZ;
/* GPU kernels start */
/// Relaksasi menggunakan penyelesaian iteratif Red-Black Gauss-Seidel (bagian Red)
///
/// \param VPotential float* Array potensial
/// \param RhoChargeDensity float* Array rapat arus
/// \param RRow int Jumlah baris di arah sumbu \f$ r \f$
/// \param ZColumn int Jumlah kolom di arah sumbu \f$ z \f$
/// \param PhiSlice int Jumlah irisan di arah sumbu \f$ \phi \f$
/// \param coef1 float* Array untuk koefisien \f$ V_{x+1,y,z} \f$
/// \param coef2 float* Array untuk koefisien \f$ V_{x-1,y,z} \f$
/// \param coef3 float* Array untuk koefisien \f$ z \f$
/// \param coef4 float* Array untuk koefisien \f$ f(r,\phi,z) \f$
__global__ void relaxationGaussSeidelRed
(
float *VPotential,
@ -52,6 +63,17 @@ __global__ void relaxationGaussSeidelRed
}
}
/// Relaksasi menggunakan penyelesaian iteratif Red-Black Gauss-Seidel (bagian Black)
///
/// \param VPotential float* Array potensial
/// \param RhoChargeDensity float* Array rapat arus
/// \param RRow int Jumlah baris di arah sumbu \f$ r \f$
/// \param ZColumn int Jumlah kolom di arah sumbu \f$ z \f$
/// \param PhiSlice int Jumlah irisan di arah sumbu \f$ \phi \f$
/// \param coef1 float* Array untuk koefisien \f$ V_{x+1,y,z} \f$
/// \param coef2 float* Array untuk koefisien \f$ V_{x-1,y,z} \f$
/// \param coef3 float* Array untuk koefisien \f$ z \f$
/// \param coef4 float* Array untuk koefisien \f$ f(r,\phi,z) \f$
__global__ void relaxationGaussSeidelBlack
(
float *VPotential,
@ -93,6 +115,21 @@ __global__ void relaxationGaussSeidelBlack
}
}
/// Menghitung residu dari hasil proses relaksasi
///
/// Rumus:
///
///
/// \param VPotential float* Array potensial
/// \param RhoChargeDensity float* Array rapat arus
/// \param DeltaResidue float* Array residu
/// \param RRow int Jumlah baris di arah sumbu \f$ r \f$
/// \param ZColumn int Jumlah kolom di arah sumbu \f$ z \f$
/// \param PhiSlice int Jumlah irisan di arah sumbu \f$ \phi \f$
/// \param coef1 float* Array untuk koefisien \f$ V_{x+1,y,z} \f$
/// \param coef2 float* Array untuk koefisien \f$ V_{x-1,y,z} \f$
/// \param coef3 float* Array untuk koefisien \f$ z \f$
/// \param icoef4 float* Array untuk koefisien invers dari \f$ f(r,\phi,z) \f$
__global__ void residueCalculation
(
float *VPotential,
@ -131,6 +168,15 @@ __global__ void residueCalculation
}
}
/// Restriksi dari finer grid ke coarser grid dengan operator Half Weighting
///
/// \f$ I_h^{2h} = \frac{1}{8} \begin{bmatrix}[ccc] 0 & 1 & 0 \\ 1 & 4 & 1\\ 0 & 1 & 0 \end{bmatrix}
///
/// \param RhoChargeDensity float* Array rapat arus
/// \param DeltaResidue float* Array residu hasil relaksasi
/// \param RRow const int Jumlah baris di arah sumbu \f$ r \f$
/// \param ZColumn const int Jumlah kolom di arah sumbu \f$ z \f$
/// \param PhiSlice const int Jumlah irisan di arah sumbu \f$ \phi \f$
__global__ void restriction2DHalf
(
float *RhoChargeDensity,
@ -168,7 +214,15 @@ __global__ void restriction2DHalf
0.125 * (DeltaResidue[finer_index_left] + DeltaResidue[finer_index_right] + DeltaResidue[finer_index_up] + DeltaResidue[finer_index_down]);
}
}
/// Restriksi dari finer grid ke coarser grid dengan operator Full Weighting
///
/// \f$ I_h^{2h} = \frac{1}{16} \begin{bmatrix}[ccc] 1 & 2 & 1 \\ 2 & 4 & 2\\ 1 & 2 & 1 \end{bmatrix}
///
/// \param RhoChargeDensity float* Array rapat arus
/// \param DeltaResidue float* Array residu hasil relaksasi
/// \param RRow const int Jumlah baris di arah sumbu \f$ r \f$
/// \param ZColumn const int Jumlah kolom di arah sumbu \f$ z \f$
/// \param PhiSlice const int Jumlah irisan di arah sumbu \f$ \phi \f$
__global__ void restriction2DFull
(
float *RhoChargeDensity,

Loading…
Cancel
Save