import qualified Algebra.ToRational as ToRational ( C )
import Data.Vector.Fixed ( Arity, S )
import NumericPrelude hiding ( abs )
-import qualified Prelude as P
+import Prelude ()
import Integration.Gaussian ( gaussian )
import Linear.Matrix (
combine ci ni = ci*(ni xi)
--- energy_true :: (Arity m, Arity n, Arity l,
--- Algebraic.C a, Eq a, RealField.C a, ToRational.C a)
--- => PDE a
--- -> Params m n (S l) a
--- -> (a -> a) -- ^ True solution @u@
--- -> (a -> a) -- ^ Derivative of true solution @u'@
--- -> a
--- energy_true pde params u u' =
--- case (bdy pde) of
--- Left (Dirichlet (x1,x2)) ->
--- sqrt $ bilinear_form u u' u u'
--- where
--- two = fromInteger 2
--- q = affine_inv (x1,x2)
--- bilinear_form w w' v v' = (x2 - x1)*(gaussian integrand)/two
--- where
--- integrand x = ((big_A pde) (q x))*(w' (q x))*(v' (q x))
--- + ((c pde) (q x))*(w (q x))*(v (q x))
-
--- _ -> error "Neumann BCs not implemented."
-
-
energy_fem :: (Arity m, Arity n, Arity l,
Algebraic.C a, Eq a, RealField.C a, ToRational.C a)
- => PDE a
- -> Params m n (S l) a
- -> a
+ => PDE a
+ -> Params m n (S l) a
+ -> a
energy_fem pde params =
(coefficients pde params) `dot` (big_F pde params)
-> a -- ^ The point @x@ at which to compute the error.
-> a
relative_error_pointwise pde params u x =
- cent * ( u_exact - u_fem ) / u_exact
+ cent * ( abs $ (u x) - u_fem ) / ( abs $ u x )
where
- u_exact = abs $ u x
u_fem = evaluate' (solution pde params) x
cent = fromInteger 100