diff --git a/Project.toml b/Project.toml
index 5b28db2..d422127 100644
--- a/Project.toml
+++ b/Project.toml
@@ -7,7 +7,7 @@ FastGaussQuadrature = "442a2c76-b920-505d-bb47-c5924d526838"
 LagrangePolynomials = "fb83b288-ce1a-45f4-939e-c34c34051e49"
 
 [sources]
-LagrangePolynomials = {rev = "main", url = "https://github.com/moment-kinetics/LagrangePolynomials"}
+LagrangePolynomials = {rev = "main", url = "https://github.com/moment-kinetics/LagrangePolynomials.jl.git"}
 
 [compat]
 LagrangePolynomials = "0.1.0"
diff --git a/README.md b/README.md
index cafa829..db3fb2b 100644
--- a/README.md
+++ b/README.md
@@ -10,13 +10,13 @@ This package computes matrices of the following form.
 A_{ij} = \int^{x_u}_{x_l} P_i(x) Q_j(x) x^p d x,
 ```
 where $`P_i`$ and $`Q_i`$ are from the set of Lagrange polynomials basis
-functions (and their first derivatives) used for  
+functions (and their first derivatives) used for
 $`C^0`$ continuous Galerkin
 finite element models, $`x_u`$ and $`x_l`$
 are the upper and lower limits of the element in the
 coordinate $`x`$, and $`p`$ is an integer.
 
-The basis functions we use here are 
+The basis functions we use here are
 ```math
 \Phi(x) = \Theta(x-x_l)\Theta(x_u - x)l_j(z(x)),
 ```
@@ -33,7 +33,7 @@ is required for the problem of interest.
 To compute the matrix $`A_{ij}`$, we use the reference coordinate
 $`z`$ to write
 ```math
-A_{ij} = \int^{1}_{-1} P_i(z) Q_j(z) (s z + c)^p d z. 
+A_{ij} = \int^{1}_{-1} P_i(z) Q_j(z) (s z + c)^p d z.
 ```
 We can also compute nonlinear matrices which have the form
 ```math
@@ -45,8 +45,8 @@ N_{ijk} = \int^{x_u}_{x_l} P_i(x) Q_j(x) R_k(x) x^p d x = \int^{1}_{-1} P_i(z) Q
 To create these matrices, the user must supply the key grid
 information for the element of interest, i.e., the set
 $`\{z_j\}`$, and the factors $`s`$ and $`c`$.
-This is done in the following function call, defining a 
-type `element_coordinates`. We use [FastGaussQuadrature](https://juliaapproximation.github.io/FastGaussQuadrature.jl/stable/)
+This is done in the following function call, defining a
+type `ElementCoordinates`. We use [FastGaussQuadrature](https://juliaapproximation.github.io/FastGaussQuadrature.jl/stable/)
 to make a set of reference nodes, but any set of nodes with good
 interpolation properties would be adequate.
 ```
@@ -57,9 +57,9 @@ x, w = gausslobatto(ngrid)
 scale = 1.0 # change these to match the size and position of elements in the problem
 shift = 0.0
 # the coordinate information
-coordinate = element_coordinates(x,scale,shift)
+coordinate = ElementCoordinates(x,scale,shift)
 ```
-Once an instance of `element_coordinates` is available, 
+Once an instance of `ElementCoordinates` is available,
 we can compute the higher-order finite element matrices
 as follows. We use the enums
 `lagrange_x`, and `d_lagrange_dx` and the integer `p` to choose the form of the desired matrix. For a mass matrix, with
@@ -67,13 +67,13 @@ Jacobian $`1`$ for example,
 we can use the following commands.
 ```
 using FiniteElementMatrices
-coordinate = element_coordinates(x,scale,shift)
-M = finite_element_matrix(lagrange_x,lagrange_x,0,coordinate)    
+coordinate = ElementCoordinates(x,scale,shift)
+M = finite_element_matrix(lagrange_x,lagrange_x,0,coordinate)
 ```
 For a stiffness matrix describing a second derivative, we could
 use the following.
 ```
-K = -finite_element_matrix(d_lagrange_dx,d_lagrange_dx,0,coordinate)    
+K = -finite_element_matrix(d_lagrange_dx,d_lagrange_dx,0,coordinate)
 ```
 The second derivative of a function $`f`$ can be found by solving the system
 ```math
@@ -86,8 +86,8 @@ in $`x`$.
 For a system with Jacobian $`x`$ (e.g., cylindrical polar coordinates), we could form the appropriate matrices to evaluate
 the 1D, radial Laplacian via
 ```
-coordinate = element_coordinates(x,scale,shift)
-M = finite_element_matrix(lagrange_x,lagrange_x,1,coordinate)    
+coordinate = ElementCoordinates(x,scale,shift)
+M = finite_element_matrix(lagrange_x,lagrange_x,1,coordinate)
 K = -finite_element_matrix(d_lagrange_dx,d_lagrange_dx,1,coordinate)
 ```
 i.e.,
@@ -96,7 +96,7 @@ i.e.,
 ```
 Finally, we can form a nonlinear stiffness matrix by inserting an extra enum argument, as follows.
 ```
-N = finite_element_matrix(lagrange_x,lagrange_x,lagrange_x,1,coordinate)    
+N = finite_element_matrix(lagrange_x,lagrange_x,lagrange_x,1,coordinate)
 ```
 
 # Examples
diff --git a/src/FiniteElementMatrices.jl b/src/FiniteElementMatrices.jl
index 9170997..e52fb5e 100644
--- a/src/FiniteElementMatrices.jl
+++ b/src/FiniteElementMatrices.jl
@@ -4,24 +4,24 @@ module FiniteElementMatrices
 
 using LagrangePolynomials: lagrange_poly,
                            lagrange_poly_derivative,
-                           lagrange_poly_data
+                           LagrangePolyData
 using FastGaussQuadrature: gausslegendre
 
 export lagrange_x,
        d_lagrange_dx,
        finite_element_matrix,
-       element_coordinates
+       ElementCoordinates
 
-@enum lagrange_function_type begin
+@enum LagrangeFunctionType begin
     lagrange_x
     d_lagrange_dx
 end
 
-struct element_coordinates
+struct ElementCoordinates
     # precomputed data for calculating
     # the Lagrange polynomials, including
     # the reference x nodes
-    lpoly_data::lagrange_poly_data
+    lpoly_data::LagrangePolyData
     # scale factor s for transform
     # v = s x + c
     # from reference grid x on [-1,1] (or (-1,1], [-1,1), (-1,1))
@@ -37,19 +37,19 @@ struct element_coordinates
     `shift` : normalisation factor, such that the physical coordinate
             `v = `scale * x_nodes + shift` in this element.
     """
-    function element_coordinates(x_nodes::AbstractArray{Float64,1},
+    function ElementCoordinates(x_nodes::AbstractArray{Float64,1},
                                     scale::Float64,
                                     shift::Float64)
         # initialise and save Lagrange polynomial
         # data for this 1D element.
-        lpoly_data = lagrange_poly_data(x_nodes)
+        lpoly_data = LagrangePolyData(x_nodes)
         return new(lpoly_data,
                     scale,
                     shift)
     end
 end
 
-function select_lagrange_function(fn_type::lagrange_function_type)
+function select_lagrange_function(fn_type::LagrangeFunctionType)
     if fn_type == lagrange_x
         func = lagrange_poly
     elseif fn_type == d_lagrange_dx
@@ -58,7 +58,7 @@ function select_lagrange_function(fn_type::lagrange_function_type)
     return func
 end
 
-function select_lagrange_prefactor(fn_type::lagrange_function_type,
+function select_lagrange_prefactor(fn_type::LagrangeFunctionType,
                                     scale::Float64)
     if fn_type == lagrange_x
         prefactor = 1.0
@@ -82,10 +82,10 @@ function get_polyz(power::Int64,
 end
 
 function finite_element_matrix(
-    fn1_type::lagrange_function_type,
-    fn2_type::lagrange_function_type,
+    fn1_type::LagrangeFunctionType,
+    fn2_type::LagrangeFunctionType,
     power::Int64,
-    coordinate::element_coordinates
+    coordinate::ElementCoordinates
     )
     lpoly_data = coordinate.lpoly_data
     ngrid = length(coordinate.lpoly_data.x_nodes)
@@ -124,11 +124,11 @@ function finite_element_matrix(
 end
 
 function finite_element_matrix(
-    fn1_type::lagrange_function_type,
-    fn2_type::lagrange_function_type,
-    fn3_type::lagrange_function_type,
+    fn1_type::LagrangeFunctionType,
+    fn2_type::LagrangeFunctionType,
+    fn3_type::LagrangeFunctionType,
     power::Int64,
-    coordinate::element_coordinates
+    coordinate::ElementCoordinates
     )
     lpoly_data = coordinate.lpoly_data
     ngrid = length(coordinate.lpoly_data.x_nodes)
diff --git a/test/runtests.jl b/test/runtests.jl
index 683c813..894d219 100644
--- a/test/runtests.jl
+++ b/test/runtests.jl
@@ -4,7 +4,7 @@ using Test: @testset, @test
 using FiniteElementMatrices: lagrange_x,
                              d_lagrange_dx,
                              finite_element_matrix,
-                             element_coordinates
+                             ElementCoordinates
 using FastGaussQuadrature: gausslobatto, gaussradau, gausslegendre
 using LinearAlgebra: lu, ldiv!, mul!
 
@@ -45,7 +45,7 @@ function test_first_derivative(;nodes::node_type=GLL,
     x = reference_nodes(nodes,ngrid)
     scale = 0.5*(y1-y0)
     shift = 0.5*(y0+y1)
-    coordinate = element_coordinates(x,scale,shift)
+    coordinate = ElementCoordinates(x,scale,shift)
     M = finite_element_matrix(lagrange_x,lagrange_x,0,coordinate)
     P = finite_element_matrix(lagrange_x,d_lagrange_dx,pp,coordinate)
     S = finite_element_matrix(d_lagrange_dx,lagrange_x,pp,coordinate)
@@ -86,7 +86,7 @@ function test_first_derivative_nonlinear_operators(
     x = reference_nodes(nodes,ngrid)
     scale = 0.5*(y1-y0)
     shift = 0.5*(y0+y1)
-    coordinate = element_coordinates(x,scale,shift)
+    coordinate = ElementCoordinates(x,scale,shift)
     M = finite_element_matrix(lagrange_x,lagrange_x,0,coordinate)
     Y100 = finite_element_matrix(d_lagrange_dx,lagrange_x,lagrange_x,pp,coordinate)
     Y010 = finite_element_matrix(lagrange_x,d_lagrange_dx,lagrange_x,pp,coordinate)
@@ -131,7 +131,7 @@ function test_second_derivative(;nodes::node_type=GLL,
     x = reference_nodes(nodes,ngrid)
     scale = 0.5*(y1-y0)
     shift = 0.5*(y0+y1)
-    coordinate = element_coordinates(x,scale,shift)
+    coordinate = ElementCoordinates(x,scale,shift)
     M = finite_element_matrix(lagrange_x,lagrange_x,0,coordinate)
     S = -finite_element_matrix(d_lagrange_dx,lagrange_x,0,coordinate)
     K = -finite_element_matrix(d_lagrange_dx,d_lagrange_dx,0,coordinate)
@@ -189,7 +189,7 @@ function test_nonlinear_operators(;nodes::node_type=GLL,
     x = reference_nodes(nodes,ngrid)
     scale = 0.5*(y1-y0)
     shift = 0.5*(y0+y1)
-    coordinate = element_coordinates(x,scale,shift)
+    coordinate = ElementCoordinates(x,scale,shift)
     M = finite_element_matrix(lagrange_x,lagrange_x,0,coordinate)
     Y000 = finite_element_matrix(lagrange_x,lagrange_x,lagrange_x,0,coordinate)
     Y100 = finite_element_matrix(d_lagrange_dx,lagrange_x,lagrange_x,0,coordinate)