Skip to content

Commit c6910b7

Browse files
pancettapancetta
committed
playing around
1 parent 493b73b commit c6910b7

File tree

8 files changed

+503
-8
lines changed

8 files changed

+503
-8
lines changed

‎pySDC/playgrounds/PETSc/playground_fft.py‎

Lines changed: 41 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,41 @@
1+
from petsc4py import PETSc
2+
3+
4+
def main():
5+
6+
n = 4
7+
da = PETSc.DMDA().create([n, n], stencil_width=1)
8+
9+
rank = PETSc.COMM_WORLD.getRank()
10+
11+
x = da.createGlobalVec()
12+
xa = da.getVecArray(x)
13+
(xs, xe), (ys, ye) = da.getRanges()
14+
print(da.getRanges())
15+
for i in range(xs, xe):
16+
for j in range(ys, ye):
17+
xa[i, j] = j*n + i
18+
print('x=', rank, x.getArray(), xs, xe, ys, ye)
19+
20+
A = da.createMatrix()
21+
A.setType(PETSc.Mat.Type.FFTW) # sparse
22+
A.setFromOptions()
23+
24+
# Istart, Iend = A.getOwnershipRange()
25+
# for I in range(Istart, Iend):
26+
# A[I, I] = 1.0
27+
28+
# communicate off-processor values
29+
# and setup internal data structures
30+
# for performing parallel operations
31+
A.assemblyBegin()
32+
A.assemblyEnd()
33+
34+
res = da.createGlobalVec()
35+
A.mult(x, res)
36+
print(rank, res.getArray())
37+
print((res-x).norm())
38+
39+
40+
if __name__ == "__main__":
41+
main()

‎pySDC/playgrounds/PETSc/playground_hypre.py‎

Lines changed: 63 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,63 @@
1+
import petsc4py, sys
2+
# args = "-ksp_type cg -pc_type hypre -pc_hypre_type boomeramg -ksp_converged_reason"
3+
# args = "-ksp_type gmres -pc_type gamg -ksp_converged_reason -ksp_monitor"
4+
args = "-ksp_type richardson -pc_type gamg -ksp_converged_reason -ksp_monitor"
5+
petsc4py.init(args)
6+
7+
from petsc4py import PETSc
8+
9+
10+
# grid size and spacing
11+
for grid_sp in [500]:
12+
m, n = grid_sp, grid_sp
13+
hx = 1.0/(m-1)
14+
hy = 1.0/(n-1)
15+
16+
# create sparse matrix
17+
A = PETSc.Mat()
18+
A.create(PETSc.COMM_WORLD)
19+
A.setSizes([m*n, m*n])
20+
A.setType('aij') # sparse
21+
A.setPreallocationNNZ(5)
22+
23+
# precompute values for setting
24+
# diagonal and non-diagonal entries
25+
diagv = 2.0/hx**2 + 2.0/hy**2
26+
offdx = -1.0/hx**2
27+
offdy = -1.0/hy**2
28+
29+
# loop over owned block of rows on this
30+
# processor and insert entry values
31+
Istart, Iend = A.getOwnershipRange()
32+
print(Istart, Iend)
33+
for I in range(Istart, Iend) :
34+
A[I,I] = diagv
35+
i = I//n # map row number to
36+
j = I - i*n # grid coordinates
37+
if i> 0 : J = I-n; A[I,J] = offdx
38+
if i< m-1: J = I+n; A[I,J] = offdx
39+
if j> 0 : J = I-1; A[I,J] = offdy
40+
if j< n-1: J = I+1; A[I,J] = offdy
41+
42+
# communicate off-processor values
43+
# and setup internal data structures
44+
# for performing parallel operations
45+
A.assemblyBegin()
46+
A.assemblyEnd()
47+
# print(A.isSymmetric())
48+
49+
# create linear solver
50+
ksp = PETSc.KSP()
51+
ksp.create(PETSc.COMM_WORLD)
52+
# obtain sol & rhs vectors
53+
x, b = A.createVecs()
54+
print(type(x))
55+
exit()
56+
x.set(0)
57+
b.set(1)
58+
# and next solve
59+
ksp.setOperators(A)
60+
ksp.setFromOptions()
61+
ksp.solve(b, x)
62+
# x.set(0)
63+
# ksp.solveTranspose(b, x)

‎pySDC/playgrounds/PETSc/playground_parallel.py‎

Lines changed: 12 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -4,13 +4,25 @@
44
def main():
55

66
n = 4
7+
8+
v = PETSc.Vec()
9+
print(type(v))
10+
v.createMPI(n, comm=PETSc.COMM_WORLD)
11+
print(type(v))
12+
# v.setSizes(n)
13+
# v.assemble()
14+
print(v.getLocalSize())
15+
exit()
16+
17+
718
da = PETSc.DMDA().create([n, n], stencil_width=1)
819

920
rank = PETSc.COMM_WORLD.getRank()
1021

1122
x = da.createGlobalVec()
1223
xa = da.getVecArray(x)
1324
(xs, xe), (ys, ye) = da.getRanges()
25+
print(da.getRanges())
1426
for i in range(xs, xe):
1527
for j in range(ys, ye):
1628
xa[i, j] = j*n + i

‎pySDC/playgrounds/VSDC/__init__.py‎

Whitespace-only changes.

‎pySDC/playgrounds/VSDC/hamiltonian_and_energy_output.py‎

Lines changed: 73 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,73 @@
1+
from pySDC.core.Hooks import hooks
2+
3+
4+
class hamiltonian_and_energy_output(hooks):
5+
6+
def __init__(self):
7+
"""
8+
Initialization of particles output
9+
"""
10+
super(hamiltonian_and_energy_output, self).__init__()
11+
self.ham_init = None
12+
self.energy_init = None
13+
14+
def pre_run(self, step, level_number):
15+
# some abbreviations
16+
L = step.levels[0]
17+
P = L.prob
18+
super(hamiltonian_and_energy_output, self).pre_run(step, level_number)
19+
self.ham_init = P.eval_hamiltonian(L.u[0])
20+
self.energy_init = P.eval_mode_energy(L.u[0])
21+
22+
def post_iteration(self, step, level_number):
23+
"""
24+
Overwrite standard post iteration hook
25+
26+
Args:
27+
step (pySDC.Step.step): the current step
28+
level_number (int): the current level number
29+
"""
30+
super(hamiltonian_and_energy_output, self).post_iteration(step, level_number)
31+
32+
# some abbreviations
33+
L = step.levels[0]
34+
P = L.prob
35+
36+
L.sweep.compute_end_point()
37+
H = P.eval_hamiltonian(L.uend)
38+
E = P.eval_mode_energy(L.uend)
39+
40+
self.add_to_stats(process=step.status.slot, time=L.time, level=-1, iter=step.status.iter,
41+
sweep=L.status.sweep, type='hamiltonian', value=H)
42+
43+
self.add_to_stats(process=step.status.slot, time=L.time, level=-1, iter=step.status.iter,
44+
sweep=L.status.sweep, type='err_hamiltonian', value=abs(self.ham_init - H))
45+
46+
self.add_to_stats(process=step.status.slot, time=L.time, level=-1, iter=step.status.iter,
47+
sweep=L.status.sweep, type='energy_iter', value=E)
48+
49+
return None
50+
51+
def post_step(self, step, level_number):
52+
"""
53+
Overwrite standard post iteration hook
54+
55+
Args:
56+
step (pySDC.Step.step): the current step
57+
level_number (int): the current level number
58+
"""
59+
super(hamiltonian_and_energy_output, self).post_step(step, level_number)
60+
61+
# some abbreviations
62+
L = step.levels[0]
63+
P = L.prob
64+
65+
E = P.eval_mode_energy(L.uend)
66+
67+
self.add_to_stats(process=step.status.slot, time=L.time, level=-1, iter=step.status.iter,
68+
sweep=L.status.sweep, type='position', value=L.uend.pos)
69+
70+
self.add_to_stats(process=step.status.slot, time=L.time, level=-1, iter=step.status.iter,
71+
sweep=L.status.sweep, type='energy_step', value=E)
72+
73+
return None

‎pySDC/playgrounds/VSDC/hamiltonian_output.py‎

Lines changed: 61 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,61 @@
1+
from pySDC.core.Hooks import hooks
2+
3+
4+
class hamiltonian_output(hooks):
5+
6+
def __init__(self):
7+
"""
8+
Initialization of particles output
9+
"""
10+
super(hamiltonian_output, self).__init__()
11+
self.ham_init = None
12+
13+
def pre_run(self, step, level_number):
14+
# some abbreviations
15+
L = step.levels[0]
16+
P = L.prob
17+
super(hamiltonian_output, self).pre_run(step, level_number)
18+
self.ham_init = P.eval_hamiltonian(L.u[0])
19+
20+
def post_iteration(self, step, level_number):
21+
"""
22+
Overwrite standard post iteration hook
23+
24+
Args:
25+
step (pySDC.Step.step): the current step
26+
level_number (int): the current level number
27+
"""
28+
super(hamiltonian_output, self).post_iteration(step, level_number)
29+
30+
# some abbreviations
31+
L = step.levels[0]
32+
P = L.prob
33+
34+
L.sweep.compute_end_point()
35+
H = P.eval_hamiltonian(L.uend)
36+
37+
self.add_to_stats(process=step.status.slot, time=L.time, level=-1, iter=step.status.iter,
38+
sweep=L.status.sweep, type='hamiltonian', value=H)
39+
40+
self.add_to_stats(process=step.status.slot, time=L.time, level=-1, iter=step.status.iter,
41+
sweep=L.status.sweep, type='err_hamiltonian', value=abs(self.ham_init - H))
42+
43+
return None
44+
45+
def post_step(self, step, level_number):
46+
"""
47+
Overwrite standard post iteration hook
48+
49+
Args:
50+
step (pySDC.Step.step): the current step
51+
level_number (int): the current level number
52+
"""
53+
super(hamiltonian_output, self).post_step(step, level_number)
54+
55+
# some abbreviations
56+
L = step.levels[0]
57+
58+
self.add_to_stats(process=step.status.slot, time=L.time, level=-1, iter=step.status.iter,
59+
sweep=L.status.sweep, type='position', value=L.uend.pos)
60+
61+
return None

0 commit comments

Comments
 (0)