3
$\begingroup$

In the finite element method, $x$-$y$ space is represented in terms of an isoparametric $\xi$-$\eta$ space using shape functions $N_{i}$ as interpolants of nodal coordinates $(x_{i},y_{i})$, where

$x=\sum_{i=1}^{n}N_{i}x_{i}$

and in this case $n=5$.

I am trying to create $x$ as a function of $\xi$ and $\eta$ and display the function as well as a simplified version of the function. See the code below.

(* Interpolants *)
N1[ξ_, η_] := -1/4 ξ (1 - ξ) (1 - η);
N2[ξ_, η_] := 1/2 (1 - ξ) (1 + ξ) (1 - η);
N3[ξ_, η_] := 1/4 ξ (1 + ξ) (1 - η);
N4[ξ_, η_] := 1/4 (1 + ξ) (1 + η);
N5[ξ_, η_] := 1/4 (1 - ξ) (1 + η);
(* Nodal Positions *)
X1 = 0; Y1 = 0;
X2 = 1/2; Y2 = 1/8;
X3 = 1; Y3 = 1/2;
X4 = 1; Y4 = 1;
X5 = 0; Y5 = 1;
x[ξ_, η_] := N1 X1 + N2 X2 + N3 X3 + N4 X4 + N5 X5;
Print[x[ξ, η]];
Print[Simplify[x[ξ, η]]];

Which produces the output

N2 / 2 + N3 + N4
N2 / 2 + N3 + N4

which is not in terms of $\xi$ and $\eta$. Once this is working I would also like to print out $y(\xi,\eta)$ which is calculated in a similar manner.

Any help would be appreciated.

$\endgroup$
7
$\begingroup$

I think the problem is how you called your N, you need to pass them the arguments as well. Like this

x[ξ_, η_] := N1[ξ, η] X1 + N2[ξ, η] X2 + N3[ξ, η] X3 +    N4[ξ, η] X4 + N5[ξ, η] X5;

And now

x[ξ, η]

Mathematica graphics

And simplifying the above gives

Mathematica graphics

Because you defined your $N_i$ as functions N1[ξ_, η_] so you need to also call them the same way you defined them.

To do what you had before, then your code should have been like this

N1 = -1/4 ξ (1 - ξ) (1 - η);
N2 = 1/2 (1 - ξ) (1 + ξ) (1 - η);
N3 = 1/4 ξ (1 + ξ) (1 - η);
N4 = 1/4 (1 + ξ) (1 + η);
N5 = 1/4 (1 - ξ) (1 + η);
(*Nodal Positions*)
X1  = 0; Y1 = 0;
X2 = 1/2; Y2 = 1/8;
X3 = 1; Y3 = 1/2;
X4 = 1; Y4 = 1;
X5 = 0; Y5 = 1;
x = N1 X1 + N2 X2 + N3 X3 + N4 X4 + N5 X5;

Mathematica graphics

| improve this answer | |
$\endgroup$
4
$\begingroup$

For your consideration:

Attributes[passdown] = {HoldFirst};

passdown[LHS : _[par : __Pattern] := RHS_] := 
 SetDelayed @@ 
  Join[Hold[LHS], 
   Hold[par][[All, 1]] /. _[p__] :> 
     Replace[Hold[RHS], s_Symbol /; DownValues[s] =!= {} :> s[p], ∞]]

Usage:

passdown[
 x[ξ_, η_] := N1 X1 + N2 X2 + N3 X3 + N4 X4 + N5 X5
]

The definition created:

?x
Global`x

x[ξ_, η_] := N1[ξ, η] X1 + N2[ξ, η] X2 + N3[ξ, η] X3 + N4[ξ, η] X4 + N5[ξ, η] X5
x[a, b]
% // Simplify
1/4 (1 - a) (1 + a) (1 - b) + 1/4 a (1 + a) (1 - b) + 1/4 (1 + a) (1 + b)

(1 + a)/2
| improve this answer | |
$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.