mulimoen
/
SummationByParts
1
0
Fork 0
Code Issues Pull Requests Packages Projects Releases Wiki Activity

zero-pad diffxi kernel

This commit is contained in:
Magnus Ulimoen 2021-03-15 19:57:50 +01:00
parent e43e71a4d8
commit 17ab18e953
4 changed files with 41 additions and 15 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

23
sbp/src/operators/algos/constmatrix.rs
View File

@ -113,7 +113,6 @@ impl<T, const N: usize> RowVector<T, N> {
}
impl<T, const M: usize, const P: usize> Matrix<T, M, P> {
#[inline(always)]
pub fn matmul_into<const N: usize>(&mut self, lhs: &Matrix<T, M, N>, rhs: &Matrix<T, N, P>)
where
T: Default + Copy + core::ops::Mul<Output = T> + core::ops::Add<Output = T>,
@ -302,6 +301,28 @@ impl<const M: usize, const N: usize> Matrix<super::Float, M, N> {
}
m
}
/// Zero extends if larger than self
pub(crate) const fn resize<const M2: usize, const N2: usize>(
&self,
) -> Matrix<super::Float, M2, N2> {
let mut m = Matrix::new([[0.0; N2]; M2]);
let m_min = if M < M2 { M } else { M2 };
let n_min = if N < N2 { N } else { N2 };
let mut i = 0;
while i < m_min {
let mut j = 0;
while j < n_min {
m.data[i][j] = self.data[i][j];
j += 1;
}
i += 1;
}
m
}
}
mod flipping {

8
sbp/src/operators/traditional4.rs
View File

@ -18,16 +18,16 @@ impl SBP4 {
1.0 / 12.0, -2.0 / 3.0, 0.0, 2.0 / 3.0, -1.0 / 12.0
]]);
#[rustfmt::skip]
const DIFF_BLOCK: Matrix<Float, 4, 6> = Matrix::new([
const MINIMAL_DIFF_BLOCK: Matrix<Float, 4, 6> = Matrix::new([
[-24.0/17.0, 59.0/34.0, -4.0/17.0, -3.0/34.0, 0.0, 0.0],
[-1.0/2.0, 0.0, 1.0/2.0, 0.0, 0.0, 0.0],
[4.0/43.0, -59.0/86.0, 0.0, 59.0/86.0, -4.0/43.0, 0.0],
[3.0/98.0, 0.0, -59.0/98.0, 0.0, 32.0/49.0, -4.0/49.0]
]);
const DIFF_BLOCKEND: super::Matrix<Float, 4, 6> =
Self::DIFF_BLOCK.flip_lr().flip_ud().flip_sign();
const DIFF_BLOCK: Matrix<Float, 4, 8> = Self::MINIMAL_DIFF_BLOCK.resize();
const DIFF_BLOCKEND: Matrix<Float, 4, 8> = Self::DIFF_BLOCK.flip_lr().flip_ud().flip_sign();
const DIFF: BlockMatrix<Float, 4, 6, 5> =
const DIFF: BlockMatrix<Float, 4, 8, 5> =
BlockMatrix::new(Self::DIFF_BLOCK, Self::DIFF_DIAG, Self::DIFF_BLOCKEND);
#[rustfmt::skip]

15
sbp/src/operators/upwind4.rs
View File

@ -14,35 +14,38 @@ impl Upwind4 {
49.0 / 144.0, 61.0 / 48.0, 41.0 / 48.0, 149.0 / 144.0
]);
#[rustfmt::skip]
const DIFF_BLOCK: Matrix<Float, 4, 7> = Matrix::new([
const MINIMAL_DIFF_BLOCK: Matrix<Float, 4, 7> = Matrix::new([
[ -72.0 / 49.0, 187.0 / 98.0, -20.0 / 49.0, -3.0 / 98.0, 0.0, 0.0, 0.0],
[-187.0 / 366.0, 0.0, 69.0 / 122.0, -16.0 / 183.0, 2.0 / 61.0, 0.0, 0.0],
[ 20.0 / 123.0, -69.0 / 82.0, 0.0, 227.0 / 246.0, -12.0 / 41.0, 2.0 / 41.0, 0.0],
[ 3.0 / 298.0, 16.0 / 149.0, -227.0 / 298.0, 0.0, 126.0 / 149.0, -36.0 / 149.0, 6.0 / 149.0],
]);
const DIFF_BLOCK: Matrix<Float, 4, 8> = Self::MINIMAL_DIFF_BLOCK.resize();
#[rustfmt::skip]
const DIFF_DIAG: RowVector<Float, 7> = RowVector::new([[
-1.0 / 24.0, 1.0 / 4.0, -7.0 / 8.0, 0.0, 7.0 / 8.0, -1.0 / 4.0, 1.0 / 24.0
]]);
const DIFF_BLOCKEND: Matrix<Float, 4, 7> = Self::DIFF_BLOCK.flip_lr().flip_ud().flip_sign();
const DIFF_BLOCKEND: Matrix<Float, 4, 8> = Self::DIFF_BLOCK.flip_lr().flip_ud().flip_sign();
const DIFF: BlockMatrix<Float, 4, 7, 7> =
const DIFF: BlockMatrix<Float, 4, 8, 7> =
BlockMatrix::new(Self::DIFF_BLOCK, Self::DIFF_DIAG, Self::DIFF_BLOCKEND);
#[rustfmt::skip]
const DISS_BLOCK: Matrix<Float, 4, 7> = Matrix::new([
const MINIMAL_DISS_BLOCK: Matrix<Float, 4, 7> = Matrix::new([
[-3.0 / 49.0, 9.0 / 49.0, -9.0 / 49.0, 3.0 / 49.0, 0.0, 0.0, 0.0],
[ 3.0 / 61.0, -11.0 / 61.0, 15.0 / 61.0, -9.0 / 61.0, 2.0 / 61.0, 0.0, 0.0],
[-3.0 / 41.0, 15.0 / 41.0, -29.0 / 41.0, 27.0 / 41.0, -12.0 / 41.0, 2.0 / 41.0, 0.0],
[3.0 / 149.0, -27.0 / 149.0, 81.0 / 149.0, -117.0 / 149.0, 90.0 / 149.0, -36.0 / 149.0, 6.0 / 149.0],
]);
const DISS_BLOCK: Matrix<Float, 4, 8> = Self::MINIMAL_DISS_BLOCK.resize();
#[rustfmt::skip]
const DISS_DIAG: RowVector<Float, 7> = Matrix::new([[
1.0 / 24.0, -1.0 / 4.0, 5.0 / 8.0, -5.0 / 6.0, 5.0 / 8.0, -1.0 / 4.0, 1.0 / 24.0
]]);
const DISS_BLOCKEND: Matrix<Float, 4, 7> = Self::DISS_BLOCK.flip_lr().flip_ud();
const DISS_BLOCKEND: Matrix<Float, 4, 8> = Self::DISS_BLOCK.flip_lr().flip_ud();
const DISS: BlockMatrix<Float, 4, 7, 7> =
const DISS: BlockMatrix<Float, 4, 8, 7> =
BlockMatrix::new(Self::DISS_BLOCK, Self::DISS_DIAG, Self::DISS_BLOCKEND);
}

10
sbp/src/operators/upwind4h2.rs
View File

@ -14,35 +14,37 @@ impl Upwind4h2 {
0.91e2/0.720e3, 0.325e3/0.384e3, 0.595e3/0.576e3, 0.1909e4/0.1920e4,
]);
#[rustfmt::skip]
const DIFF_BLOCK: Matrix<Float, 4, 7> = Matrix::new([
const MINIMAL_DIFF_BLOCK: Matrix<Float, 4, 7> = Matrix::new([
[-3.95604395604396e+00, 5.41758241758242e+00, -1.94505494505495e+00, 4.83516483516484e-01, 0.0, 0.0, 0.0],
[-8.09025641025641e-01, 0.00000000000000e+00, 1.03948717948718e+00, -2.47384615384615e-01, 1.69230769230769e-02, 0.0, 0.0],
[2.37983193277311e-01, -8.51680672268908e-01, 0.00000000000000e+00, 7.35966386554622e-01, -1.36134453781513e-01, 1.38655462184874e-02, 0.0],
[-6.14632442814737e-02, 2.10581456259822e-01, -7.64623712240265e-01, 0.00000000000000e+00, 7.42535358826611e-01, -1.41435306443164e-01, 1.44054478784704e-02],
]);
const DIFF_BLOCK: Matrix<Float, 4, 8> = Self::MINIMAL_DIFF_BLOCK.resize();
#[rustfmt::skip]
const DIFF_DIAG: RowVector<Float, 7> = RowVector::new([[
-1.43229166666667e-02, 1.40625000000000e-01, -7.38281250000000e-01, 0.00000000000000e+00, 7.38281250000000e-01, -1.40625000000000e-01, 1.43229166666667e-02
]]);
const DIFF: BlockMatrix<Float, 4, 7, 7> = BlockMatrix::new(
const DIFF: BlockMatrix<Float, 4, 8, 7> = BlockMatrix::new(
Self::DIFF_BLOCK,
Self::DIFF_DIAG,
Self::DIFF_BLOCK.flip_lr().flip_ud().flip_sign(),
);
#[rustfmt::skip]
const DISS_BLOCK: Matrix<Float, 4, 7> = Matrix::new([
const MINIMAL_DISS_BLOCK: Matrix<Float, 4, 7> = Matrix::new([
[-2.76989342315976e-01, 5.19355016842454e-01, -3.46236677894969e-01, 1.03871003368491e-01, 0.0, 0.0, 0.0],
[7.75570158484731e-02, -1.62342481638964e-01, 1.47715500579822e-01, -7.98531117124082e-02, 1.69230769230769e-02, 0.0, 0.0],
[-4.23630758836198e-02, 1.21027405937249e-01, -1.91609307039399e-01, 1.82272708078206e-01, -8.31932773109244e-02, 1.38655462184874e-02, 0.0],
[1.32037874021759e-02, -6.79734450144910e-02, 1.89370108794365e-01, -2.78654929966754e-01, 2.16081718177056e-01, -8.64326872708224e-02, 1.44054478784704e-02],
]);
const DISS_BLOCK: Matrix<Float, 4, 8> = Self::MINIMAL_DISS_BLOCK.resize();
#[rustfmt::skip]
const DISS_DIAG: RowVector<Float, 7> = RowVector::new([[
1.43229166666667e-02, -8.59375000000000e-02, 2.14843750000000e-01, -2.86458333333333e-01, 2.14843750000000e-01, -8.59375000000000e-02, 1.43229166666667e-02,
]]);
const DISS: BlockMatrix<Float, 4, 7, 7> = BlockMatrix::new(
const DISS: BlockMatrix<Float, 4, 8, 7> = BlockMatrix::new(
Self::DISS_BLOCK,
Self::DISS_DIAG,
Self::DISS_BLOCK.flip_lr().flip_ud(),