Page Menu
Home
c4science
Search
Configure Global Search
Log In
Files
F88607709
anisotropic_square.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
Sat, Oct 19, 17:39
Size
4 KB
Mime Type
text/x-python
Expires
Mon, Oct 21, 17:39 (1 d, 23 h)
Engine
blob
Format
Raw Data
Handle
21797288
Attached To
R6746 RationalROMPy
anisotropic_square.py
View Options
import
numpy
as
np
import
matplotlib.pyplot
as
plt
from
itertools
import
product
from
anisotropic_square_engine
import
(
AnisotropicSquareEngine
as
engine
,
AnisotropicSquareEnginePoles
as
plsEx
)
from
rrompy.sampling
import
SamplingEngineStandard
as
SES
from
rrompy.reduction_methods
import
(
NearestNeighbor
as
NN
,
RationalInterpolantGreedyPivotedGreedy
as
RIGPG
)
from
rrompy.parameter.parameter_sampling
import
QuadratureSampler
as
QS
from
rrompy.parameter
import
localSparseGrid
as
LSG
zs
,
Ls
=
[
10.
,
50.
],
[
.
2
,
1.2
]
z0
,
L0
,
n
=
np
.
mean
(
zs
),
np
.
mean
(
Ls
),
50
murange
=
[[
zs
[
0
],
Ls
[
0
]],
[
zs
[
-
1
],
Ls
[
-
1
]]]
np
.
random
.
seed
(
4020
)
mu
=
[
zs
[
0
]
+
np
.
random
.
rand
()
*
(
zs
[
-
1
]
-
zs
[
0
]),
Ls
[
0
]
+
np
.
random
.
rand
()
*
(
Ls
[
-
1
]
-
Ls
[
0
])]
solver
=
engine
(
z0
,
L0
,
n
)
fighandles
=
[]
params
=
{
"POD"
:
True
,
"nTestPoints"
:
100
,
"greedyTol"
:
1e-4
,
"S"
:
3
,
"polybasisMarginal"
:
"MONOMIAL_WENDLAND"
,
"polybasis"
:
"LEGENDRE"
,
'samplerPivot'
:
QS
(
zs
,
"UNIFORM"
),
'trainSetGenerator'
:
QS
(
zs
,
"UNIFORM"
),
'errorEstimatorKind'
:
"LOOK_AHEAD_RES"
,
"MMarginal"
:
2
,
"SMarginal"
:
3
,
"greedyTolMarginal"
:
1e-2
,
"samplerMarginalGrid"
:
LSG
(
Ls
),
"radialDirectionalWeightsMarginal"
:
[
2.
],
"matchingWeight"
:
1.
}
for
tol
,
shared
in
product
([
1.
,
3.
],
[
1.
,
0.
]):
print
(
"Testing cutoff tolerance {} with shared ratio {}."
.
format
(
tol
,
shared
))
params
[
'cutOffTolerance'
]
=
tol
params
[
'cutOffSharedRatio'
]
=
shared
approx
=
RIGPG
([
0
],
solver
,
mu0
=
[
z0
,
L0
],
approx_state
=
True
,
approxParameters
=
params
,
verbosity
=
5
)
approx
.
setupApprox
(
"LAST"
)
print
(
"--- Approximant ---"
)
approx
.
plotApprox
(
mu
,
name
=
'u_app'
)
approx
.
plotHF
(
mu
,
name
=
'u_HF'
)
approx
.
plotErr
(
mu
,
name
=
'err_app'
)
approx
.
plotRes
(
mu
,
name
=
'res_app'
)
normErr
=
approx
.
normErr
(
mu
)[
0
]
normSol
=
approx
.
normHF
(
mu
)[
0
]
normRes
=
approx
.
normRes
(
mu
)[
0
]
normRHS
=
approx
.
normRHS
(
mu
)[
0
]
print
(
"SolNorm:
\t
{:.5e}
\n
Err_app:
\t
{:.5e}
\n
ErrRel_app:
\t
{:.5e}"
.
format
(
normSol
,
normErr
,
normErr
/
normSol
))
print
(
"RHSNorm:
\t
{:.5e}
\n
Res_app:
\t
{:.5e}
\n
ResRel_app:
\t
{:.5e}"
.
format
(
normRHS
,
normRes
,
normRes
/
normRHS
))
print
(
"--- Closest snapshot ---"
)
eng
=
SES
(
solver
,
verbosity
=
0
)
eng
.
nsamples
=
approx
.
samplingEngine
.
nsamplesCoalesced
eng
.
mus
=
approx
.
samplingEngine
.
musCoalesced
eng
.
samples
=
approx
.
samplingEngine
.
samples_fullCoalesced
paramsNN
=
{
'S'
:
eng
.
nsamples
,
'POD'
:
False
,
'sampler'
:
QS
(
murange
,
"UNIFORM"
)}
approxNN
=
NN
(
solver
,
mu0
=
[
z0
,
L0
],
approx_state
=
True
,
approxParameters
=
paramsNN
,
verbosity
=
0
)
approxNN
.
setSamples
(
eng
)
approxNN
.
plotApprox
(
mu
,
name
=
'u_close'
)
approxNN
.
plotHF
(
mu
,
name
=
'u_HF'
)
approxNN
.
plotErr
(
mu
,
name
=
'err_close'
)
approxNN
.
plotRes
(
mu
,
name
=
'res_close'
)
normErr
=
approxNN
.
normErr
(
mu
)[
0
]
normSol
=
approxNN
.
normHF
(
mu
)[
0
]
normRes
=
approxNN
.
normRes
(
mu
)[
0
]
normRHS
=
approxNN
.
normRHS
(
mu
)[
0
]
print
(
"SolNorm:
\t
{:.5e}
\n
Err_close:
\t
{:.5e}
\n
ErrRel_close:
\t
{:.5e}"
.
format
(
normSol
,
normErr
,
normErr
/
normSol
))
print
(
"RHSNorm:
\t
{:.5e}
\n
Res_close:
\t
{:.5e}
\n
ResRel_close:
\t
{:.5e}"
.
format
(
normRHS
,
normRes
,
normRes
/
normRHS
))
verb
=
approx
.
verbosity
approx
.
verbosity
=
0
tspace
=
np
.
linspace
(
Ls
[
0
],
Ls
[
-
1
],
100
)
for
j
,
t
in
enumerate
(
tspace
):
plsE
=
plsEx
(
t
,
0.
,
zs
[
-
1
])
pls
=
approx
.
getPoles
([
None
,
t
])
pls
[
np
.
abs
(
np
.
imag
(
pls
))
>
1e-5
]
=
np
.
nan
if
j
==
0
:
polesE
=
np
.
empty
((
len
(
tspace
),
len
(
plsE
)))
poles
=
np
.
empty
((
len
(
tspace
),
len
(
pls
)))
polesE
[:]
=
np
.
nan
if
len
(
plsE
)
>
polesE
.
shape
[
1
]:
nanR
=
np
.
empty
((
len
(
tspace
),
len
(
plsE
)
-
polesE
.
shape
[
1
]))
nanR
[:]
=
np
.
nan
polesE
=
np
.
hstack
((
polesE
,
nanR
))
polesE
[
j
,
:
len
(
plsE
)]
=
np
.
real
(
plsE
)
poles
[
j
]
=
np
.
real
(
pls
)
approx
.
verbosity
=
verb
fighandles
+=
[
plt
.
figure
(
figsize
=
(
17
,
5
))]
ax1
=
fighandles
[
-
1
]
.
add_subplot
(
1
,
2
,
1
)
ax2
=
fighandles
[
-
1
]
.
add_subplot
(
1
,
2
,
2
)
ax1
.
plot
(
poles
,
tspace
)
ax1
.
set_ylim
(
Ls
)
ax1
.
set_xlabel
(
'mu_1'
)
ax1
.
set_ylabel
(
'mu_2'
)
ax1
.
grid
()
ax2
.
plot
(
polesE
,
tspace
,
'k-.'
,
linewidth
=
1
)
ax2
.
plot
(
poles
,
tspace
)
for
mm
in
approx
.
musMarginal
:
ax2
.
plot
(
zs
,
[
mm
[
0
,
0
]]
*
2
,
'k--'
,
linewidth
=
1
)
ax2
.
set_xlim
(
zs
)
ax2
.
set_ylim
(
Ls
)
ax2
.
set_xlabel
(
'mu_1'
)
ax2
.
set_ylabel
(
'mu_2'
)
ax2
.
grid
()
plt
.
show
()
print
(
"
\n
"
)
Event Timeline
Log In to Comment