Broken again.

src/JuliaFEM.jl

@@ -7,7 +7,10 @@ using Lexicon
 using Logging
 @Logging.configure(level=DEBUG)
 
+Logging.info("loading types")
 include("types.jl") # type definitions
+Logging.info("loading elements")
+include("elements.jl") # elements
 include("math.jl") # basic mathematical operations
 
 include("elasticity_solver.jl")
@@ -15,4 +18,6 @@ include("xdmf.jl")
 include("abaqus_reader.jl")
 include("interfaces.jl")
 
+export set_coordinates, get_coordinates, set_material
+
 end # module

src/elements.jl

@@ -0,0 +1,92 @@
+# This file is a part of JuliaFEM.
+# License is MIT: see https://github.com/JuliaFEM/JuliaFEM.jl/blob/master/LICENSE.md
+
+abstract CG <: Element # Lagrange element family
+abstract Quad4 <: CG # 4 node quadrangle elements
+
+
+"""
+Evaluate basis functions in point xi.
+"""
+function get_basis(el::Quad4)
+    (xi) -> [(1-xi[1])*(1-xi[2])/4
+             (1+xi[1])*(1-xi[2])/4
+             (1+xi[1])*(1+xi[2])/4
+             (1-xi[1])*(1+xi[2])/4]
+end
+function get_basis(el::Quad4, xi)
+    get_basis(el)(xi)
+end
+
+"""
+Evaluate partial derivatives of basis function w.r.t
+dimensionless coordinate xi, i.e. dbasis/dxi
+"""
+function get_dbasisdxi(el::Quad4)
+    (xi) -> [-(1-xi[2])/4.0    -(1-xi[1])/4.0
+              (1-xi[2])/4.0    -(1+xi[1])/4.0
+              (1+xi[2])/4.0     (1+xi[1])/4.0
+             -(1+xi[2])/4.0     (1-xi[1])/4.0]
+end
+function get_dbasisdxi(el::Quad4, xi)
+    get_dbasisdxi(el)(xi)
+end
+
+
+
+"""
+Get jacobian of element evaluated at point xi
+"""
+function get_jacobian(el::Element, xi)
+    dbasisdxi = get_dbasisdxi(el)
+    X = get_coordinates(el)
+    J = interpolate(X, dbasisdxi, xi)'
+    return J
+end
+
+"""
+Evaluate partial derivatives of basis function w.r.t
+material description X, i.e. dbasis/dX
+"""
+function get_dbasisdX(el::CG)
+    function get_dbasisdX_(xi)
+        dbasisdxi = get_dbasisdxi(el, xi)
+        J = get_jacobian(el, xi)
+        dbasisdxi*inv(J)
+    end
+end
+function get_dbasisdX(el::CG, xi)
+    get_dbasisdX(el)(xi)
+end
+
+"""
+Return coordinates of element in array of size dim x nnodes
+"""
+function get_coordinates(el::Element)
+    # Make sure you define at least this field to your element if you want
+    # to build everything yourself
+    el.coordinates
+end
+
+function set_coordinates(el::Element, coordinates)
+    el.coordinates = coordinates
+end
+
+function set_material(el::Element, lambda, mu)
+    el.attributes["lambda"] = lambda
+    el.attributes["mu"] = mu
+end
+
+"""
+Get element id
+"""
+function get_element_id(el::Element)
+    el.id
+end
+
+function get_integration_points(el::Element)
+    el.integration_points
+end
+
+
+

src/math.jl

@@ -60,26 +60,11 @@ function interpolate{T<:Real}(field::Array{T,2}, basis::Function, ip)
     return result
 end
 function interpolate(e::Element, field::ASCIIString, x::Array{Float64,1}; derivative=false)
-    return interpolate(e.attributes[field], derivative ? e.shape_functions.dbasis : e.shape_functions.basis, x)
+    basis = derivative ? get_dbasisdxi(e) : get_basis(e)
+    return interpolate(e.attributes[field], basis, x)
 end
 
 
-"""
-
-"""
-function get_basis(el::Element, xi)
-    return el.shape_functions.basis(xi)
-end
-
-"""
-Return partial derivatives of shape functions w.r.t X using chain rule.
-"""
-function get_dbasisdX(el::Element, xi)
-    J = interpolate(el, "coordinates", xi; derivative=true)
-    dbasisdX = el.shape_functions.dbasis(xi)*inv(J')
-    return dbasisdX
-end
-
 
 """
 Linearize function f w.r.t some given field, i.e. calculate dR/du