#! /usr/bin/env python # def fem1d_bvp_linear_test00 ( ): #*****************************************************************************80 # ## FEM1D_BVP_LINEAR_TEST00 tests FEM1D_BVP_LINEAR. # # Discussion: # # - uxx + u = x for 0 < x < 1 # u(0) = u(1) = 0 # # exact = x - sinh(x) / sinh(1) # exact' = 1 - cosh(x) / sinh(1) # # Location: # # http://people.sc.fsu.edu/~jburkardt/m_src/fem1d_bvp_linear/fem1d_bvp_linear_test00.py # # Licensing: # # This code is distributed under the GNU LGPL license. # # Modified: # # 09 January 2015 # # Author: # # John Burkardt # # Reference: # # Dianne O'Leary, # Scientific Computing with Case Studies, # SIAM, 2008, # ISBN13: 978-0-898716-66-5, # LC: QA401.O44. # import matplotlib.pyplot as plt import numpy as np import platform from fem1d_bvp_linear import fem1d_bvp_linear from h1s_error_linear import h1s_error_linear from l1_error import l1_error from l2_error_linear import l2_error_linear from max_error_linear import max_error_linear n = 11 print ( '' ) print ( 'FEM1D_BVP_LINEAR_TEST00' ) print ( ' Python version: %s' % ( platform.python_version ( ) ) ) print ( ' Solve -( A(x) U\'(x) )\' + C(x) U(x) = F(x)' ) print ( ' for 0 < x < 1, with U(0) = U(1) = 0.' ) print ( ' A(X) = 1.0' ) print ( ' C(X) = 1.0' ) print ( ' F(X) = X' ) print ( ' U(X) = X - SINH(X) / SINH(1)' ) print ( '' ) print ( ' Number of nodes = %d' % ( n ) ) # # Geometry definitions. # x_lo = 0.0 x_hi = 1.0 x = np.linspace ( x_lo, x_hi, n ) u = fem1d_bvp_linear ( n, a00, c00, f00, x ) g = np.zeros ( n ) for i in range ( 0, n ): g[i] = exact00 ( x[i] ) # # Print a table. # print ( '' ) print ( ' I X U Uexact Error' ) print ( '' ) for i in range ( 0, n ): print ( ' %4d %8f %8f %8f %8e' \ % ( i, x[i], u[i], g[i], abs ( u[i] - g[i] ) ) ) # # Compute error norms. # e1 = l1_error ( n, x, u, exact00 ) e2 = l2_error_linear ( n, x, u, exact00 ) h1s = h1s_error_linear ( n, x, u, exactp00 ) mx = max_error_linear ( n, x, u, exact00 ) print ( '' ) print ( ' l1 norm of error = %g' % ( e1 ) ) print ( ' L2 norm of error = %g' % ( e2 ) ) print ( ' Seminorm of error = %g' % ( h1s ) ) print ( ' Max norm of error = %g' % ( mx ) ) # # Plot the computed solution. # fig = plt.figure ( ) plt.plot ( x, u, 'bo-' ) plt.xlabel ( '<---X--->' ) plt.ylabel ( '<---U(X)--->' ) plt.title ( 'FEM1D_BVP_LINEAR_TEST00 Solution' ) plt.savefig ( 'fem1d_bvp_linear_test00.png' ) # plt.show ( ) # # Terminate. # print ( '' ) print ( 'FEM1D_BVP_LINEAR_TEST00' ) print ( ' Normal end of execution.' ) return def a00 ( x ): value = 1.0 return value def c00 ( x ): value = 1.0 return value def exact00 ( x ): from math import sinh value = x - sinh ( x ) / sinh ( 1.0 ) return value def exactp00 ( x ): from math import cosh from math import sinh value = 1.0 - cosh ( x ) / sinh ( 1.0 ) return value def f00 ( x ): value = x return value if ( __name__ == '__main__' ): from timestamp import timestamp timestamp ( ) fem1d_bvp_linear_test00 ( ) timestamp ( )