|
106 | 106 | .. autoclass:: LineOfCompressionKernel |
107 | 107 | :show-inheritance: |
108 | 108 | :members: mapper_method |
| 109 | +.. autoclass:: HeatKernel |
| 110 | + :show-inheritance: |
| 111 | + :members: mapper_method |
109 | 112 |
|
110 | 113 | Derivatives |
111 | 114 | ----------- |
@@ -1119,6 +1122,76 @@ def get_pde_as_diff_op(self) -> LinearPDESystemOperator: |
1119 | 1122 | return laplacian(w) |
1120 | 1123 |
|
1121 | 1124 |
|
| 1125 | +class HeatKernel(ExpressionKernel): |
| 1126 | + r"""The Green's function for the heat equation given by |
| 1127 | + :math:`e^{-r^2/{4 \alpha t}}/\sqrt{(4 \pi \alpha t)^d}` |
| 1128 | + where :math:`d` is the number of spatial dimensions. |
| 1129 | +
|
| 1130 | + .. note:: |
| 1131 | +
|
| 1132 | + This kernel cannot be used in an FMM yet and can only |
| 1133 | + be used in expansions and evaluations that occur forward |
| 1134 | + in the time dimension. |
| 1135 | + """ |
| 1136 | + heat_alpha_name: str |
| 1137 | + |
| 1138 | + mapper_method: ClassVar[str] = "map_heat_kernel" |
| 1139 | + |
| 1140 | + def __init__(self, spatial_dims: int, heat_alpha_name: str = "alpha"): |
| 1141 | + dim = spatial_dims + 1 |
| 1142 | + d = make_sym_vector("d", dim) |
| 1143 | + t = d[-1] |
| 1144 | + r = pymbolic_real_norm_2(d[:-1]) |
| 1145 | + alpha = SpatialConstant(heat_alpha_name) |
| 1146 | + expr = var("exp")(-r**2/(4 * alpha * t)) / var("sqrt")(t**(dim - 1)) |
| 1147 | + scaling = 1/var("sqrt")((4*var("pi")*alpha)**(dim - 1)) |
| 1148 | + |
| 1149 | + super().__init__( |
| 1150 | + dim, |
| 1151 | + expression=expr, |
| 1152 | + global_scaling_const=scaling, |
| 1153 | + ) |
| 1154 | + |
| 1155 | + self.heat_alpha_name = heat_alpha_name |
| 1156 | + |
| 1157 | + @property |
| 1158 | + @override |
| 1159 | + def is_complex_valued(self) -> bool: |
| 1160 | + return False |
| 1161 | + |
| 1162 | + def update_persistent_hash(self, key_hash, key_builder): |
| 1163 | + key_hash.update(type(self).__name__.encode("utf8")) |
| 1164 | + key_builder.rec(key_hash, (self.dim, self.heat_alpha_name)) |
| 1165 | + |
| 1166 | + @override |
| 1167 | + def __repr__(self): |
| 1168 | + return f"HeatKnl{self.dim - 1}D" |
| 1169 | + |
| 1170 | + @override |
| 1171 | + def get_args(self): |
| 1172 | + return [ |
| 1173 | + KernelArgument( |
| 1174 | + loopy_arg=lp.ValueArg(self.heat_alpha_name, np.float64), |
| 1175 | + )] |
| 1176 | + |
| 1177 | + def get_derivative_taker(self, dvec, rscale, sac): |
| 1178 | + """Return a :class:`sumpy.derivative_taker.ExprDerivativeTaker` instance |
| 1179 | + that supports taking derivatives of the base kernel with respect to dvec. |
| 1180 | + """ |
| 1181 | + from sumpy.derivative_taker import ExprDerivativeTaker |
| 1182 | + return ExprDerivativeTaker(self.get_expression(dvec), dvec, rscale, |
| 1183 | + sac) |
| 1184 | + |
| 1185 | + @override |
| 1186 | + def get_pde_as_diff_op(self): |
| 1187 | + from sumpy.expansion.diff_op import diff, laplacian, make_identity_diff_op |
| 1188 | + alpha = sym.Symbol(self.heat_alpha_name) |
| 1189 | + w = make_identity_diff_op(self.dim - 1, time_dependent=True) |
| 1190 | + time_diff = [0]*self.dim |
| 1191 | + time_diff[-1] = 1 |
| 1192 | + return diff(w, tuple(time_diff)) - laplacian(w) * alpha |
| 1193 | + |
| 1194 | + |
1122 | 1195 | # }}} |
1123 | 1196 |
|
1124 | 1197 |
|
@@ -1590,6 +1663,7 @@ def map_expression_kernel(self, kernel: ExpressionKernel) -> Kernel: |
1590 | 1663 | map_elasticity_kernel = map_expression_kernel |
1591 | 1664 | map_line_of_compression_kernel = map_expression_kernel |
1592 | 1665 | map_stresslet_kernel = map_expression_kernel |
| 1666 | + map_heat_kernel = map_expression_kernel |
1593 | 1667 |
|
1594 | 1668 | def map_axis_target_derivative(self, kernel: AxisTargetDerivative) -> Kernel: |
1595 | 1669 | return type(kernel)(kernel.axis, self.rec(kernel.inner_kernel)) |
@@ -1662,6 +1736,7 @@ def map_expression_kernel(self, kernel: ExpressionKernel) -> int: |
1662 | 1736 | map_yukawa_kernel = map_expression_kernel |
1663 | 1737 | map_line_of_compression_kernel = map_expression_kernel |
1664 | 1738 | map_stresslet_kernel = map_expression_kernel |
| 1739 | + map_heat_kernel = map_expression_kernel |
1665 | 1740 |
|
1666 | 1741 | @override |
1667 | 1742 | def map_axis_target_derivative(self, kernel: AxisTargetDerivative) -> int: |
|
0 commit comments