Page Menu
Home
c4science
Search
Configure Global Search
Log In
Files
F60186079
vector_problem_engine_base.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
Sun, Apr 28, 03:43
Size
8 KB
Mime Type
text/x-python
Expires
Tue, Apr 30, 03:43 (1 d, 23 h)
Engine
blob
Format
Raw Data
Handle
17314467
Attached To
R6746 RationalROMPy
vector_problem_engine_base.py
View Options
# Copyright (C) 2018 by the RROMPy authors
#
# This file is part of RROMPy.
#
# RROMPy is free software: you can redistribute it and/or modify
# it under the terms of the GNU Lesser General Public License as published by
# the Free Software Foundation, either version 3 of the License, or
# (at your option) any later version.
#
# RROMPy is distributed in the hope that it will be useful,
# but WITHOUT ANY WARRANTY; without even the implied warranty of
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
# GNU Lesser General Public License for more details.
#
# You should have received a copy of the GNU Lesser General Public License
# along with RROMPy. If not, see <http://www.gnu.org/licenses/>.
#
import
fenics
as
fen
import
numpy
as
np
from
matplotlib
import
pyplot
as
plt
from
rrompy.utilities.base.types
import
Np1D
,
strLst
from
rrompy.utilities.base
import
purgeList
,
getNewFilename
from
.problem_engine_base
import
ProblemEngineBase
__all__
=
[
'VectorProblemEngineBase'
]
class
VectorProblemEngineBase
(
ProblemEngineBase
):
"""
Generic solver for parametric vector problems.
Attributes:
verbosity: Verbosity level.
BCManager: Boundary condition manager.
V: Real FE space.
u: Generic trial functions for variational form evaluation.
v: Generic test functions for variational form evaluation.
As: Scipy sparse array representation (in CSC format) of As.
bs: Numpy array representation of bs.
energyNormMatrix: Scipy sparse matrix representing inner product over
V.
bsmu: Mu value of last bs evaluation.
liftDirichletDatamu: Mu value of last Dirichlet datum evaluation.
liftedDirichletDatum: Dofs of Dirichlet datum lifting.
mu0BC: Mu value of last Dirichlet datum lifting.
degree_threshold: Threshold for ufl expression interpolation degree.
"""
nAs
,
nbs
=
1
,
1
def
__init__
(
self
,
degree_threshold
:
int
=
np
.
inf
,
verbosity
:
int
=
10
,
timestamp
:
bool
=
True
):
super
()
.
__init__
(
degree_threshold
=
degree_threshold
,
verbosity
=
verbosity
,
timestamp
=
timestamp
)
self
.
V
=
fen
.
VectorFunctionSpace
(
fen
.
UnitSquareMesh
(
10
,
10
),
"P"
,
1
)
def
plot
(
self
,
u
:
Np1D
,
name
:
str
=
"u"
,
save
:
str
=
None
,
what
:
strLst
=
'all'
,
saveFormat
:
str
=
"eps"
,
saveDPI
:
int
=
100
,
show
:
bool
=
True
,
**
figspecs
):
"""
Do some nice plots of the complex-valued function with given dofs.
Args:
u: numpy complex array with function dofs.
name(optional): Name to be shown as title of the plots. Defaults to
'u'.
save(optional): Where to save plot(s). Defaults to None, i.e. no
saving.
what(optional): Which plots to do. If list, can contain 'ABS',
'PHASE', 'REAL', 'IMAG'. If str, same plus wildcard 'ALL'.
Defaults to 'ALL'.
saveFormat(optional): Format for saved plot(s). Defaults to "eps".
saveDPI(optional): DPI for saved plot(s). Defaults to.
show(optional): Whether to show figure. Defaults to True.
figspecs(optional key args): Optional arguments for matplotlib
figure creation.
"""
if
isinstance
(
what
,
(
str
,)):
if
what
.
upper
()
==
'ALL'
:
what
=
[
'ABS'
,
'PHASE'
,
'REAL'
,
'IMAG'
]
else
:
what
=
[
what
]
what
=
purgeList
(
what
,
[
'ABS'
,
'PHASE'
,
'REAL'
,
'IMAG'
],
listname
=
self
.
name
()
+
".what"
,
baselevel
=
1
)
if
'figsize'
not
in
figspecs
.
keys
():
figspecs
[
'figsize'
]
=
(
13.
*
max
(
len
(
what
),
1
)
/
4
,
3
)
if
len
(
what
)
>
0
:
for
j
in
range
(
self
.
V
.
num_sub_spaces
()):
subplotcode
=
100
+
len
(
what
)
*
10
II
=
self
.
V
.
sub
(
j
)
.
dofmap
()
.
dofs
()
Vj
=
self
.
V
.
sub
(
j
)
.
collapse
()
plt
.
figure
(
**
figspecs
)
plt
.
jet
()
if
'ABS'
in
what
:
uAb
=
fen
.
Function
(
Vj
)
uAb
.
vector
()
.
set_local
(
np
.
abs
(
u
[
II
]))
subplotcode
=
subplotcode
+
1
plt
.
subplot
(
subplotcode
)
p
=
fen
.
plot
(
uAb
,
title
=
"|{}_comp{}|"
.
format
(
name
,
j
))
plt
.
colorbar
(
p
)
if
'PHASE'
in
what
:
uPh
=
fen
.
Function
(
Vj
)
uPh
.
vector
()
.
set_local
(
np
.
angle
(
u
[
II
]))
subplotcode
=
subplotcode
+
1
plt
.
subplot
(
subplotcode
)
p
=
fen
.
plot
(
uPh
,
title
=
"phase({}_comp{})"
.
format
(
name
,
j
))
plt
.
colorbar
(
p
)
if
'REAL'
in
what
:
uRe
=
fen
.
Function
(
Vj
)
uRe
.
vector
()
.
set_local
(
np
.
real
(
u
[
II
]))
subplotcode
=
subplotcode
+
1
plt
.
subplot
(
subplotcode
)
p
=
fen
.
plot
(
uRe
,
title
=
"Re({}_comp{})"
.
format
(
name
,
j
))
plt
.
colorbar
(
p
)
if
'IMAG'
in
what
:
uIm
=
fen
.
Function
(
Vj
)
uIm
.
vector
()
.
set_local
(
np
.
imag
(
u
[
II
]))
subplotcode
=
subplotcode
+
1
plt
.
subplot
(
subplotcode
)
p
=
fen
.
plot
(
uIm
,
title
=
"Im({}_comp{})"
.
format
(
name
,
j
))
plt
.
colorbar
(
p
)
if
save
is
not
None
:
save
=
save
.
strip
()
plt
.
savefig
(
getNewFilename
(
"{}_comp{}_fig_"
.
format
(
save
,
j
),
saveFormat
),
format
=
saveFormat
,
dpi
=
saveDPI
)
if
show
:
plt
.
show
()
plt
.
close
()
try
:
if
len
(
what
)
>
1
:
figspecs
[
'figsize'
]
=
(
2.
/
len
(
what
)
*
figspecs
[
'figsize'
][
0
],
figspecs
[
'figsize'
][
1
])
elif
len
(
what
)
==
0
:
figspecs
[
'figsize'
]
=
(
2.
*
figspecs
[
'figsize'
][
0
],
figspecs
[
'figsize'
][
1
])
if
len
(
what
)
==
0
or
'ABS'
in
what
or
'REAL'
in
what
:
uVRe
=
fen
.
Function
(
self
.
V
)
uVRe
.
vector
()
.
set_local
(
np
.
real
(
u
))
plt
.
figure
(
**
figspecs
)
plt
.
jet
()
p
=
fen
.
plot
(
uVRe
,
title
=
"{}_Re"
.
format
(
name
),
mode
=
"displacement"
)
plt
.
colorbar
(
p
)
if
save
is
not
None
:
save
=
save
.
strip
()
plt
.
savefig
(
getNewFilename
(
"{}_disp_Re_fig_"
.
format
(
save
),
saveFormat
),
format
=
saveFormat
,
dpi
=
saveDPI
)
plt
.
show
()
plt
.
close
()
if
'ABS'
in
what
or
'IMAG'
in
what
:
uVIm
=
fen
.
Function
(
self
.
V
)
uVIm
.
vector
()
.
set_local
(
np
.
imag
(
u
))
plt
.
figure
(
**
figspecs
)
plt
.
jet
()
p
=
fen
.
plot
(
uVIm
,
title
=
"{}_Im"
.
format
(
name
),
mode
=
"displacement"
)
plt
.
colorbar
(
p
)
if
save
is
not
None
:
save
=
save
.
strip
()
plt
.
savefig
(
getNewFilename
(
"{}_disp_Im_fig_"
.
format
(
save
,
j
),
saveFormat
),
format
=
saveFormat
,
dpi
=
saveDPI
)
if
show
:
plt
.
show
()
plt
.
close
()
except
:
pass
def
plotmesh
(
self
,
name
:
str
=
"Mesh"
,
save
:
str
=
None
,
saveFormat
:
str
=
"eps"
,
saveDPI
:
int
=
100
,
show
:
bool
=
True
,
**
figspecs
):
"""
Do a nice plot of the mesh.
Args:
u: numpy complex array with function dofs.
name(optional): Name to be shown as title of the plots. Defaults to
'u'.
save(optional): Where to save plot(s). Defaults to None, i.e. no
saving.
saveFormat(optional): Format for saved plot(s). Defaults to "eps".
saveDPI(optional): DPI for saved plot(s). Defaults to 100.
show(optional): Whether to show figure. Defaults to True.
figspecs(optional key args): Optional arguments for matplotlib
figure creation.
"""
plt
.
figure
(
**
figspecs
)
fen
.
plot
(
self
.
V
.
mesh
())
if
save
is
not
None
:
save
=
save
.
strip
()
plt
.
savefig
(
getNewFilename
(
"{}_msh_"
.
format
(
save
),
saveFormat
),
format
=
saveFormat
,
dpi
=
saveDPI
)
if
show
:
plt
.
show
()
plt
.
close
()
Event Timeline
Log In to Comment