Page Menu
Home
c4science
Search
Configure Global Search
Log In
Files
F86893292
active_remeshing.py
No One
Temporary
Actions
Download File
Edit File
Delete File
View Transforms
Subscribe
Mute Notifications
Award Token
Subscribers
None
File Metadata
Details
File Info
Storage
Attached
Created
Wed, Oct 9, 06:16
Size
4 KB
Mime Type
text/x-python
Expires
Fri, Oct 11, 06:16 (2 d)
Engine
blob
Format
Raw Data
Handle
21501430
Attached To
R6746 RationalROMPy
active_remeshing.py
View Options
import
numpy
as
np
from
pickle
import
load
from
matplotlib
import
pyplot
as
plt
from
active_remeshing_engine
import
ActiveRemeshingEngine
from
rrompy.reduction_methods
import
(
RationalInterpolantGreedyPivotedGreedyNoMatch
as
RIGPGNM
,
RationalInterpolantGreedyPivotedGreedy
as
RIGPG
)
from
rrompy.parameter.parameter_sampling
import
(
QuadratureSampler
as
QS
,
SparseGridSampler
as
SGS
)
zs
,
ts
=
[
0.
,
100.
],
[
0.
,
.
5
]
z0
,
t0
,
n
=
np
.
mean
(
zs
),
np
.
mean
(
ts
),
150
solver
=
ActiveRemeshingEngine
(
z0
,
t0
,
n
)
mus
=
[[
z0
,
ts
[
0
]],
[
z0
,
ts
[
1
]]]
for
mu
in
mus
:
u
=
solver
.
solve
(
mu
,
return_state
=
True
)[
0
]
Y
=
solver
.
applyC
(
u
,
mu
)[
0
][
0
]
_
=
solver
.
plot
(
u
,
what
=
"REAL"
,
name
=
"u(z={}, t={})"
.
format
(
*
mu
),
is_state
=
True
,
figsize
=
(
12
,
4
))
print
(
"Y(z={}, t={}) = {} (solution norm squared)"
.
format
(
*
mu
,
Y
))
fighandles
=
[]
with
open
(
"./active_remeshing_hf_samples.pkl"
,
"rb"
)
as
f
:
zspace
,
tspace
,
Yex
=
load
(
f
)
# zspace = np.linspace(zs[0], zs[-1], 198)
# tspace = np.linspace(ts[0], ts[-1], 50)
# Yex = np.log10([[solver.solve([z, t]) for t in tspace] for z in zspace])
# (from a ~2h45m simulation on one node of the EPFL Helvetios cluster)
for
match
in
[
0
,
1
]:
params
=
{
"POD"
:
True
,
"S"
:
5
,
"SMarginal"
:
3
,
"greedyTol"
:
1e-4
,
"nTestPoints"
:
500
,
"polybasis"
:
"LEGENDRE"
,
"maxIterMarginal"
:
25
,
"trainSetGenerator"
:
QS
(
zs
,
"UNIFORM"
),
"samplerPivot"
:
QS
(
zs
,
"CHEBYSHEV"
),
"samplerMarginal"
:
SGS
(
ts
),
"errorEstimatorKind"
:
"LOOK_AHEAD_OUTPUT"
,
"errorEstimatorKindMarginal"
:
"LOOK_AHEAD_RECOVER"
}
if
match
:
print
(
"
\n
Testing output-based matching with weight 1."
)
params
[
"greedyTolMarginal"
]
=
1e-2
params
[
"matchingWeight"
]
=
1.
params
[
"sharedRatio"
]
=
.
75
params
[
"polybasisMarginal"
]
=
"PIECEWISE_LINEAR_UNIFORM"
algo
=
RIGPG
else
:
print
(
"
\n
Testing matching-free approach."
)
params
[
"greedyTolMarginal"
]
=
2e-2
algo
=
RIGPGNM
approx
=
algo
([
0
],
solver
,
mu0
=
[
z0
,
t0
],
approxParameters
=
params
,
verbosity
=
15
)
approx
.
setupApprox
(
"ALL"
)
verb
,
verbTM
=
approx
.
verbosity
,
approx
.
trainedModel
.
verbosity
approx
.
verbosity
,
approx
.
trainedModel
.
verbosity
=
0
,
0
for
j
,
t
in
enumerate
(
tspace
):
out
=
approx
.
getApprox
(
np
.
pad
(
zspace
.
reshape
(
-
1
,
1
),
[(
0
,
0
),
(
0
,
1
)],
"constant"
,
constant_values
=
t
))
pls
=
approx
.
getPoles
([
None
,
t
])
pls
[
np
.
abs
(
np
.
imag
(
pls
))
>
1e-5
]
=
np
.
nan
if
j
==
0
:
Ys
=
np
.
empty
((
len
(
zspace
),
len
(
tspace
)))
poles
=
np
.
empty
((
len
(
tspace
),
len
(
pls
)))
Ys
[:,
j
]
=
np
.
log10
(
out
.
data
)
if
len
(
pls
)
>
poles
.
shape
[
1
]:
poles
=
np
.
pad
(
poles
,
[(
0
,
0
),
(
0
,
len
(
pls
)
-
poles
.
shape
[
1
])],
"constant"
,
constant_values
=
np
.
nan
)
poles
[
j
,
:
len
(
pls
)]
=
np
.
real
(
pls
)
approx
.
verbosity
,
approx
.
trainedModel
.
verbosity
=
verb
,
verbTM
Ymin
,
Ymax
=
min
(
np
.
min
(
Ys
),
np
.
min
(
Yex
)),
max
(
np
.
max
(
Ys
),
np
.
max
(
Yex
))
fighandles
+=
[
plt
.
figure
(
figsize
=
(
15
,
5
))]
ax1
=
fighandles
[
-
1
]
.
add_subplot
(
1
,
2
,
1
)
ax2
=
fighandles
[
-
1
]
.
add_subplot
(
1
,
2
,
2
)
if
match
:
ax1
.
plot
(
poles
,
tspace
)
else
:
ax1
.
plot
(
poles
,
tspace
,
'k.'
)
ax1
.
set_ylim
(
ts
)
ax1
.
set_xlabel
(
'z'
)
ax1
.
set_ylabel
(
't'
)
ax1
.
grid
()
if
match
:
ax2
.
plot
(
poles
,
tspace
)
else
:
ax2
.
plot
(
poles
,
tspace
,
'k.'
)
for
mm
in
approx
.
musMarginal
:
ax2
.
plot
(
zs
,
[
mm
[
0
,
0
]]
*
2
,
'k--'
,
linewidth
=
1
)
ax2
.
set_xlim
(
zs
)
ax2
.
set_ylim
(
ts
)
ax2
.
set_xlabel
(
'z'
)
ax2
.
set_ylabel
(
't'
)
ax2
.
grid
()
plt
.
show
()
print
(
"Approximate poles"
)
fighandles
+=
[
plt
.
figure
(
figsize
=
(
15
,
5
))]
ax1
=
fighandles
[
-
1
]
.
add_subplot
(
1
,
2
,
1
)
ax2
=
fighandles
[
-
1
]
.
add_subplot
(
1
,
2
,
2
)
p
=
ax1
.
contourf
(
np
.
repeat
(
zspace
.
reshape
(
-
1
,
1
),
len
(
tspace
),
axis
=
1
),
np
.
repeat
(
tspace
.
reshape
(
1
,
-
1
),
len
(
zspace
),
axis
=
0
),
Ys
,
vmin
=
Ymin
,
vmax
=
Ymax
,
cmap
=
"gray_r"
,
levels
=
np
.
linspace
(
Ymin
,
Ymax
,
50
))
plt
.
colorbar
(
p
,
ax
=
ax1
)
ax1
.
set_xlabel
(
'z'
)
ax1
.
set_ylabel
(
't'
)
ax1
.
grid
()
p
=
ax2
.
contourf
(
np
.
repeat
(
zspace
.
reshape
(
-
1
,
1
),
len
(
tspace
),
axis
=
1
),
np
.
repeat
(
tspace
.
reshape
(
1
,
-
1
),
len
(
zspace
),
axis
=
0
),
Yex
,
vmin
=
Ymin
,
vmax
=
Ymax
,
cmap
=
"gray_r"
,
levels
=
np
.
linspace
(
Ymin
,
Ymax
,
50
))
ax2
.
set_xlabel
(
'z'
)
ax2
.
set_ylabel
(
't'
)
ax2
.
grid
()
plt
.
colorbar
(
p
,
ax
=
ax2
)
plt
.
show
()
print
(
"Approximate and exact output
\n
"
)
Event Timeline
Log In to Comment