diff --git a/README.md b/README.md index 733348c..a209e0d 100644 --- a/README.md +++ b/README.md @@ -22,7 +22,7 @@ code = create_basis( (0.5, 0.0, 0.5), # N9 (0.0, 0.5, 0.5), # N10 ), - "1 + u + v + w + u*v + v*w + w*u + u^2 + v^2 + w^2", + :(1 + u + v + w + u*v + v*w + w*u + u^2 + v^2 + w^2), ) ``` @@ -114,11 +114,11 @@ code = create_basis( ( 0.0, 0.0, 1.0), # N5 ), ( - "1/8 * (1-u) * (1-v) * (1-w)", - "1/8 * (1+u) * (1-v) * (1-w)", - "1/8 * (1+u) * (1+v) * (1-w)", - "1/8 * (1-u) * (1+v) * (1-w)", - "1/2 * (1+w)", + :(1/8 * (1-u) * (1-v) * (1-w)), + :(1/8 * (1+u) * (1-v) * (1-w)), + :(1/8 * (1+u) * (1+v) * (1-w)), + :(1/8 * (1-u) * (1+v) * (1-w)), + :(1/2 * (1+w)), ), ) eval(code) diff --git a/src/create_basis.jl b/src/create_basis.jl index 895bc5b..245c0a1 100644 --- a/src/create_basis.jl +++ b/src/create_basis.jl @@ -85,7 +85,7 @@ end """ Calculate derivatives of basis functions with respect to parameters u, v, w. """ -function calculate_interpolation_polynomial_derivatives(basis::Vector, dim::Int) +function calculate_interpolation_polynomial_derivatives(basis, dim) equations = Vector[] vars = [:u, :v, :w] nbasis = length(basis) @@ -101,7 +101,7 @@ function calculate_interpolation_polynomial_derivatives(basis::Vector, dim::Int) return equations end -function create_basis{nbasis,dim}(name::Symbol, description::String, X::NTuple{nbasis, NTuple{dim, Float64}}, basis::Vector, dbasis::Vector) +function create_basis{nbasis,dim}(name::Symbol, description::String, X::NTuple{nbasis, NTuple{dim, Float64}}, basis, dbasis) Q = Expr(:block) for i=1:nbasis @@ -177,3 +177,8 @@ function create_basis{nbasis,dim}(name::Symbol, description::String, X::NTuple{n dbasis = calculate_interpolation_polynomial_derivatives(basis, dim) return create_basis(name, description, X, basis, dbasis) end + +function create_basis{nbasis,dim}(name::Symbol, description::String, X::NTuple{nbasis, NTuple{dim, Float64}}, basis::NTuple{nbasis, Expr}) + dbasis = calculate_interpolation_polynomial_derivatives(basis, dim) + return create_basis(name, description, X, basis, dbasis) +end diff --git a/test/test_create_basis.jl b/test/test_create_basis.jl index 0060738..2936bdb 100644 --- a/test/test_create_basis.jl +++ b/test/test_create_basis.jl @@ -6,6 +6,8 @@ using FEMBasis #using FEMBasis: calculate_basis_coefficients, calculate_interpolation_polynomials, # calculate_interpolation_polynomial_derivatives, create_lagrange_basis +@testset "FEMBasis.jl" begin + @testset "Calculate interpolation polynomial matrix" begin p = :(1 + u + v) X = ((0.0, 0.0), (1.0, 0.0), (0.0, 1.0)) @@ -18,9 +20,6 @@ end p = :(1 + u + v) A = [1.0 0.0 0.0; 1.0 1.0 0.0; 1.0 0.0 1.0] equations = FEMBasis.calculate_interpolation_polynomials(p, A) - for j=1:3 - println(equations[j]) - end @test equations[1] == :(1.0 + -u + -v) @test equations[2] == :(+u) @test equations[3] == :(+v) @@ -28,7 +27,6 @@ end @testset "Calculate derivatives of interpolation polynomials" begin - println("call derivatve function") basis = [:(1.0 - u - v), :(1.0u), :(1.0v)] dbasis = FEMBasis.calculate_interpolation_polynomial_derivatives(basis, 2) @test isapprox(dbasis[1][1], -1.0) @@ -44,7 +42,6 @@ end basis = [:(1.0 - u - v), :(1.0u), :(1.0v)] dbasis = Vector[[-1.0, -1.0], [1.0, 0.0], [0.0, 1.0]] code = FEMBasis.create_basis(:TestTriangle, "test triangle", X, basis, dbasis) - println(code) eval(code) N = zeros(1,size(TestTriangle, 2)) X = get_reference_element_coordinates(TestTriangle) @@ -63,12 +60,10 @@ end @test isapprox(dN, [-1.0 1.0 0.0; -1.0 0.0 1.0]) end -@testset "test create basis, ansatz polynomial given" begin +@testset "test create basis, N given" begin X = ((0.0, 0.0), (1.0, 0.0), (0.0, 1.0)) - # p = :(1 + u + v) - p = "1 + u + v" - code = FEMBasis.create_basis(:TestTriangle, "test triangle", X, p) - println(code) + basis = ( :(1 - u - v), :(1.0u), :(1.0v) ) + code = FEMBasis.create_basis(:TestTriangle, "test triangle", X, basis) eval(code) N = zeros(1,size(TestTriangle, 2)) X = get_reference_element_coordinates(TestTriangle) @@ -86,3 +81,30 @@ end eval_dbasis!(TestTriangle, dN, X[3]) @test isapprox(dN, [-1.0 1.0 0.0; -1.0 0.0 1.0]) end + +@testset "test create basis, ansatz polynomial given" begin + X = ((0.0, 0.0), (1.0, 0.0), (0.0, 1.0)) + p1 = :(1 + u + v) + p2 = "1 + u + v" + code1 = FEMBasis.create_basis(:TestTriangle, "test triangle", X, p1) + code2 = FEMBasis.create_basis(:TestTriangle, "test triangle", X, p2) + @test code1 == code2 + eval(code1) + N = zeros(1,size(TestTriangle, 2)) + X = get_reference_element_coordinates(TestTriangle) + for i=1:length(TestTriangle) + eval_basis!(TestTriangle, N, X[i]) + N_expected = zeros(1, length(TestTriangle)) + N_expected[i] = 1.0 + @test isapprox(N, N_expected) + end + dN = zeros(size(TestTriangle)...) + eval_dbasis!(TestTriangle, dN, X[1]) + @test isapprox(dN, [-1.0 1.0 0.0; -1.0 0.0 1.0]) + eval_dbasis!(TestTriangle, dN, X[2]) + @test isapprox(dN, [-1.0 1.0 0.0; -1.0 0.0 1.0]) + eval_dbasis!(TestTriangle, dN, X[3]) + @test isapprox(dN, [-1.0 1.0 0.0; -1.0 0.0 1.0]) +end + +end