#!/usr/bin/env python
# -*- coding: utf-8 -*-
# Python GetFEM interface
#
# Copyright (C) 2020-2020 Tetsuo Koyama.
#
# This file is a part of GetFEM
#
# GetFEM  is  free software;  you  can  redistribute  it  and/or modify it
# under  the  terms  of the  GNU  Lesser General Public License as published
# by  the  Free Software Foundation;  either version 2.1 of the License,  or
# (at your option) any later version.
# This program  is  distributed  in  the  hope  that it will be useful,  but
# WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
# or  FITNESS  FOR  A PARTICULAR PURPOSE.  See the GNU Lesser General Public
# License for more details.
# You  should  have received a copy of the GNU Lesser General Public License
# along  with  this program;  if not, write to the Free Software Foundation,
# Inc., 51 Franklin St, Fifth Floor, Boston, MA  02110-1301, USA.
#
############################################################################
""" 1D Truss problem test.

  This program is used to check that python-getfem is working. This is
  also a good example of use of GetFEM.

  $Id$
"""
import getfem as gf
import numpy as np

E = 1.0E+05   # Young Modulus (N/mm^2)
A = 1.0       # Cross-Sectional Area (mm^2)
L = 1000.0    # Length (mm)

mesh = gf.Mesh("cartesian", [0, L])
mfu = gf.MeshFem(mesh, 1)

element_degree = 1
mfu.set_classical_fem(element_degree)

mim = gf.MeshIm(mesh, element_degree*2)

model = gf.Model("real")
model.add_fem_variable("u", mfu)

model.add_initialized_data("D", [E*A])
model.add_bilaplacian_brick(mim, "u", "D")

model.assembly()
SM = model.tangent_matrix()

print("Computed matrix is...")
print(SM)
print("Expected matrix is...")
SM = gf.Spmat("empty", 2, 2)
SM.add(0, 0, E*A/L)
SM.add(0, 1, -E*A/L)
SM.add(1, 0, -E*A/L)
SM.add(1, 1, E*A/L)
print(SM)