Page Menu
Home
c4science
Search
Configure Global Search
Log In
Files
F85268602
double_slit_engine.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
Fri, Sep 27, 22:03
Size
2 KB
Mime Type
text/x-python
Expires
Sun, Sep 29, 22:03 (2 d)
Engine
blob
Format
Raw Data
Handle
21126529
Attached To
R6746 RationalROMPy
double_slit_engine.py
View Options
import
numpy
as
np
import
ufl
import
fenics
as
fen
import
mshr
from
rrompy.utilities.base.decorators
import
(
affine_construct
,
nonaffine_construct
)
from
rrompy.hfengines.linear_problem
import
ScatteringProblemEngine
from
rrompy.utilities.numerical.hash_derivative
import
(
hashDerivativeToIdx
as
hashD
)
from
rrompy.solver.fenics
import
fenZERO
,
fenics2Vector
class
DoubleSlitEngine
(
ScatteringProblemEngine
):
def
__init__
(
self
,
k0
:
float
,
n
:
int
):
super
()
.
__init__
(
mu0
=
[
k0
])
self
.
_affinePoly
=
False
delta
,
eps
=
.
1
,
.
01
mesh
=
mshr
.
generate_mesh
(
mshr
.
Circle
(
fen
.
Point
(
0.
,
0.
),
5.
)
-
mshr
.
Rectangle
(
fen
.
Point
(
-
5.
,
-
delta
),
fen
.
Point
(
-.
75
,
delta
))
-
mshr
.
Rectangle
(
fen
.
Point
(
-.
5
,
-
delta
),
fen
.
Point
(
.
5
,
delta
))
-
mshr
.
Rectangle
(
fen
.
Point
(
.
75
,
-
delta
),
fen
.
Point
(
5.
,
delta
)),
n
)
self
.
V
=
fen
.
FunctionSpace
(
mesh
,
"P"
,
1
)
self
.
DirichletBoundary
=
lambda
x
,
on_boundary
:
(
on_boundary
and
np
.
abs
(
x
[
1
])
<=
delta
and
np
.
abs
(
x
[
1
])
>
delta
-
eps
)
self
.
NeumannBoundary
=
lambda
x
,
on_boundary
:
(
on_boundary
and
np
.
abs
(
x
[
1
])
<=
delta
-
eps
)
self
.
RobinBoundary
=
"REST"
def
getDirichletValues
(
self
,
mu
=
[]):
mu
=
self
.
checkParameter
(
mu
)
x
,
y
=
fen
.
SpatialCoordinate
(
self
.
V
.
mesh
())[:]
c
,
s
=
.
5
,
-
.
5
*
3.
**
.
5
muR
,
muI
=
np
.
real
(
mu
[
0
])[
0
],
np
.
imag
(
mu
[
0
])[
0
]
mod
=
-
muI
*
(
c
*
x
+
s
*
y
)
angle
=
muR
*
(
c
*
x
+
s
*
y
)
DR
=
fen
.
exp
(
mod
)
*
fen
.
cos
(
angle
)
DI
=
fen
.
exp
(
mod
)
*
fen
.
sin
(
angle
)
DR
=
ufl
.
conditional
(
ufl
.
ge
(
y
,
0
),
DR
,
fenZERO
)
DI
=
ufl
.
conditional
(
ufl
.
ge
(
y
,
0
),
DI
,
fenZERO
)
return
DR
,
DI
@affine_construct
def
A
(
self
,
mu
=
[],
der
=
0
):
return
ScatteringProblemEngine
.
A
(
self
,
mu
,
der
)
@nonaffine_construct
def
b
(
self
,
mu
=
[],
der
=
0
):
derI
=
hashD
(
der
)
if
hasattr
(
der
,
"__len__"
)
else
der
if
derI
<
0
:
return
self
.
baselineb
()
if
derI
>
0
:
raise
Exception
(
"Derivatives not implemented."
)
fen0
=
fen
.
inner
(
fenZERO
,
self
.
v
)
*
fen
.
dx
DR
,
DI
=
self
.
getDirichletValues
(
mu
)
DBCR
=
fen
.
DirichletBC
(
self
.
V
,
DR
,
self
.
DirichletBoundary
)
DBCI
=
fen
.
DirichletBC
(
self
.
V
,
DI
,
self
.
DirichletBoundary
)
return
(
fenics2Vector
(
fen0
,
{},
DBCR
,
1
)
+
1.j
*
fenics2Vector
(
fen0
,
{},
DBCI
,
1
))
Event Timeline
Log In to Comment