例題1-5 多項式 Lagrange內插法 f (x) = log (x) 求x=1.5 , 2.5 , 3.5 f(x)之值 = ? 並計算誤差值=?

例題1-5 多項式 Lagrange內插法
f (x) = ln (x)

已知　平面上4點

   x                  y=f(x)
===================
 1.0                 0.0
 2.0                 0.693
 3.0                 1.099
 4.0                  1.386

求x=1.5 , 2.5 , 3.5   f(x)之值 = ?
並計算誤差值=?



#/********** Lagrange's interpolation ***************/
import math

def f(x):
    return math.log(x)

n=4
x= [0.0 for i in range(n+1)]     #x [n] 矩陣
y= [0.0 for i in range(n+1)]     #f [n] 矩陣
x=[1.0 , 2.0  ,3.0 , 4.0]
y=[0.0 , 0.693 , 1.099 ,1.386 ]

n1=3
P=[1.5 , 2.5 , 3.5]
Px=[0.0 , 0.0 , 0.0]

#f(x)=log(x)
a=1.5   #x=1.5 求 f(1.5)= ?

print("\nEnter the number of the terms of the table: ",n)
print("\nEnter the respective values of the variables x and y: \n")
print("The values of x =")
for i  in range (0,n) :     #
    print( round( x[i],4),"\t",end='')
print()
print("The values of y =")
for i  in range (0,n) :     #
    print( round( y[i],4),"\t",end='')

print("\nEnter the value of the x to find the respective value of y ")
for i  in range (0,n1) :     #
    print( round( P[i],4),"\t",end='')
print()

for m in range (0 , n1):
    a=P[m]
    k=0.0
    for i in range(0, n):
        s=1.0
        t=1.0
        for j in range(0,n):
            if(j!=i):
                s=s*(a-x[j]);
                t=t*(x[i]-x[j]);
        k=k+((s/t)*y[i]);
    Px[m]=k
    print("\nThe {%.2f} to find respective value of the variable y is: {%.6f}" %(a, k))
 
print(" \nx\t\tP(x)\t\t\tf(x)\t\t\t|f(x)-P(x)|")
for i in range (0 , n1):
    print ("{%2.3f}\t\t{%5.6f}\t\t{%5.6f}\t\t{%5.6f}" %(P[i],Px[i], f(P[i]), abs(Px[i]-f(P[i]) ) ) )




輸出畫面
========= RESTART: F:/2018-09勤益科大數值分析/數值分析/PYTHON/EX1-5.py ============

Enter the number of the terms of the table:  4

Enter the respective values of the variables x and y: 

The values of x =
1.0 2.0 3.0 4.0
The values of y =
0.0 0.693 1.099 1.386
Enter the value of the x to find the respective value of y 
1.5 2.5 3.5

The {1.50} to find respective value of the variable y is: {0.392875}

The {2.50} to find respective value of the variable y is: {0.921375}

The {3.50} to find respective value of the variable y is: {1.246875}

x P(x) f(x) |f(x)-P(x)|
{1.500} {0.392875} {0.405465} {0.012590}
{2.500} {0.921375} {0.916291} {0.005084}
{3.500} {1.246875} {1.252763} {0.005888}
>>> 

例題1-2 Lagrange內插法 已知 平面上3點 求x=1.5 , f(x)= ?

例題1-2 Lagrange內插法

已知　平面上3點

   x                  y=f(x)
===================
 1.0                 0.0
 2.0                 0.693
 3.0                 1.099

求x=1.5 ,  f(x)= ?

#/********** Lagrange's interpolation ***************/
n=3
x= [0.0 for i in range(n+1)]     #x [n] 矩陣
y= [0.0 for i in range(n+1)]     #f [n] 矩陣
x=[1.0 , 2.0  ,3.0]
y=[0.0 , 0.693 , 1.099 ]

a=1.5   #x=1.5 求 f(1.5)= ?

print("\nEnter the number of the terms of the table: ",n)
print("\nEnter the respective values of the variables x and y: \n")
print("The values of x =")
for i  in range (0,n) :     # 
    print( round( x[i],4),"\t",end='')
print()
print("The values of y =")
for i  in range (0,n) :     # 
    print( round( y[i],4),"\t",end='')

print("\n\nEnter the value of the x to find the respective value of y ",a)

k=0.0
for i in range(0, n):
    s=1.0
    t=1.0
    for j in range(0,n):
        if(j!=i):
            s=s*(a-x[j]);
            t=t*(x[i]-x[j]);
    k=k+((s/t)*y[i]);

print("\nThe respective value of the variable y is: {%f}" %(k))




輸出畫面
======== RESTART: F:/2018-09勤益科大數值分析/數值分析/PYTHON/EX1-2.py ============

Enter the number of the terms of the table:  3

Enter the respective values of the variables x and y: 

The values of x =
1.0 2.0 3.0
The values of y =
0.0 0.693 1.099

Enter the value of the x to find the respective value of y  1.5

The respective value of the variable y is: {0.382375}
>>> 

例題1-1 Lagrange內插法 已知 平面上二點 求x=1.5 , f(x)= ?

例題1-1 Lagrange內插法
已知　平面上二點

   x                  y=f(x)
===================
 1.0                 0.0
 2.0                 0.693

求x=1.5 ,  f(x)= ?


#/********** Lagrange's interpolation ***************/
n=2
x= [0.0 for i in range(n+1)]     #x [n] 矩陣
y= [0.0 for i in range(n+1)]     #f [n] 矩陣
x=[1.0 , 2.0]
y=[0.0 , 0.693]
a=1.5   #x=1.5 求 f(1.5)= ?

k=0.0
d=1

print("\nEnter the number of the terms of the table: ",n)
print("\nEnter the respective values of the variables x and y: \n")
print("The values of x =")
for i  in range (0,n) :     #
    print( round( x[i],4),"\t",end='')
print()
print("The values of y =")
for i  in range (0,n) :     #
    print( round( y[i],4),"\t",end='')

print("\n Enter the value of the x to find the respective value of y ",a)

while(d==1):
    for i in range(0, n):
        s=1.0
        t=1.0
        for j in range(0,n):
            if(j!=i):
                s=s*(a-x[j]);
                t=t*(x[i]-x[j]);
       
        k=k+((s/t)*y[i]);

    print("\nThe respective value of the variable y is: {%f}" %(k))
    d=0

輸出畫面

========= RESTART: F:/2018-09勤益科大數值分析/數值分析/PYTHON/EX1-1.py ===========

Enter the number of the terms of the table:  2

Enter the respective values of the variables x and y: 

The values of x =

1.0 2.0

The values of y =

0.0 0.693

 Enter the value of the x to find the respective value of y  1.5

The respective value of the variable y is: {0.346500}

>>> 

範例7-5 請用有限差法 (Finite difference Method) 解 非線性常微分方程式 y'' + y' ^ 2 + y = ln(x)

範例7-5 請用有限差法 (Finite difference Method) 解
非線性常微分方程式
y'' +  y' ^ 2  + y = ln(x)

1 <= x <= 2 , y(1)= 0.0  ,  y(2)= ln(2)=0.6931472

其真實解 W(x) = ln(x)
取 h=0.1  , n=99 , h=0.01 , n=99

/* ex7-5.c uses finite difference method to solve
 * nonlinear ordinary differential equation with boundary
 * conditions, y"+p(x)y'+q(x)y=r(x),a<=x<=b,
 * y(a)=alfa, y(b)=bata. After transfer ordinary
 * differential equation into system of linear algebra
 * equations, then call function tridg() to solve
 * tridiagonal equations.
 */
#include <stdio.h>
#include <math.h>
#define  p(x,y)  (y)
#define  q(x,y)  (1.0)
#define  r(x,y)  (log(x))
#define  w(x)    (log(x))
void tridg(int,double [],double [],double [],double []);
void main()
{
   int i,k,n;
   double a[100],b[100],c[100],d[100],y[100],dy[100],
  h,x1,x,xn,aa,bb,alfa,bata,err;
   scanf("n=%d aa=%lf bb=%lf alfa=%lf bata=%lf",
  &n,&aa,&bb,&alfa,&bata);
   h=(bb-aa)/(n+1);
   for(i=1;i<=n;i++)
      y[i]=0.0;
   
   for(k=1;k<=100;k++)
   {
        x1=aa+h;
        /* dy[1]=y1'=(y2-y0)/(2h) */
        dy[1]=(1.0/(2*h))*(y[2]-alfa);
        b[1]=pow(h,2)*q((x1),(y[1]))-2.0;
        c[1]=(1+(h/2.0)*p((x1),(dy[1])));
        d[1]=pow(h,2)*r((x1),(y[1]))-(1.0-(h/2.0)*
    p((x1),(dy[1])))*alfa;
        for(i=2;i<=n-1;i++)
        {
            x=aa+i*h;
        /* dy[i]=yi' */
        dy[i]=(1.0/(2*h))*(y[i+1]-y[i-1]);
            a[i]=1-(h/2.0)*p((x),(dy[i]));
            b[i]=pow(h,2)*q((x),(y[i]))-2.0;
            c[i]=1+(h/2.0)*p((x),(dy[i]));
            d[i]=pow(h,2)*r((x),(y[i]));
        }
        xn=aa+n*h;
        /* dy[n]=yn' */
        dy[n]=(1.0/(2*h))*(bata-y[n-1]);
        a[n]=1-(h/2.0)*p((xn),(dy[n]));
        b[n]=pow(h,2)*q((xn),(y[n]))-2.0;
        d[n]=pow(h,2)*r((xn),(y[n]))-(1+(h/2.0)*
        p((xn),(dy[n])))*bata;
        tridg(n,a,b,c,d);
        err=0.0;
        for(i=1;i<=n;i++)
        err=err+fabs(d[i]-y[i]);
        if(err >0.001)
        {
        for(i=1;i<=n;i++)
            y[i]=d[i];
        }
        else
        goto bound;
   }
   bound:
   printf("The iterations=%d\n",k);
   printf("x       y(x)       w(x)      |y(x)-w(x)|\n");
   printf("%5.3lf %10.7lf %10.7lf %10.7lf\n",
   aa,alfa,w(aa),fabs(alfa-w(aa)));
   for(i=1;i<=n;i++)
   {
      x=aa+i*h;
      if(i%10==0)
      printf("%5.3lf %10.7lf %10.7lf %10.7lf\n",
      x,d[i],w(x),fabs(d[i]-w(x)));
   }
   printf("%5.3lf %10.7lf %10.7lf %10.7lf\n",
    bb,bata,w(bb),fabs(bata-w(bb)));
   return;
}
void tridg(int n,double a[],double b[],double c[],double d[])
{
   int i;
   double r;
   for(i=2;i<=n;i++)
   {
      r=a[i]/b[i-1];
      b[i]=b[i]-r*c[i-1];
      d[i]=d[i]-r*d[i-1];
   }
   /* The answers are stored in d[i] */
   d[n]=d[n]/b[n];
   for(i=n-1;i>=1;i--)
      d[i]=(d[i]-c[i]*d[i+1])/b[i];
   return;
}


輸入資料
n=99 aa=1.0 bb=2.0 alfa=0.0 bata=0.6931472

輸出資料
The iterations=6
x       y(x)       w(x)      |y(x)-w(x)|
1.000  0.0000000  0.0000000  0.0000000
1.100  0.0953101  0.0953102  0.0000001
1.200  0.1823215  0.1823216  0.0000001
1.300  0.2623644  0.2623643  0.0000001
1.400  0.3364727  0.3364722  0.0000004
1.500  0.4054658  0.4054651  0.0000007
1.600  0.4700045  0.4700036  0.0000009
1.700  0.5306291  0.5306283  0.0000009
1.800  0.5877874  0.5877867  0.0000007
1.900  0.6418543  0.6418539  0.0000004
2.000  0.6931472  0.6931472  0.0000000

範例7-4 請用有限差法 (Finite difference Method) 解 非線性常微分方程式的解 y'' + y y' - y^2 =0

範例7-4 請用有限差法 (Finite difference Method) 解
非線性常微分方程式的解
y'' + y y' - y^2 =0

1 <= x <= 2 , y(1)= 0.5  ,  y(2)= 1/3

其真實解 W(x) = 1/ (1+x) 
取 h=0.1  , n=99 , h=0.01 , n=99 


'''
/* ex7-4.py uses finite difference method to solve
 * ordinary differential equation with boundary
 * conditions, y"+p(x)y'+q(x)y=r(x),a<=x<=b,
 * y(a)=alfa, y(b)=bata. After transfer ordinary
 * differential equation into system of linear algebra
 * equations, then call function tridg() to solve
 * tridiagonal equations.
 */
 '''
import math

def  p(x,y):
    return (y)
   
def q(x,y):
    return (0.0)

def r(x,y):
    return (math.pow(y,3))

def w(x):
    return (1.0 / (1+x) )


def tridg( n , a, b,  c,  d ):
    for i in range (2 , n+1):
        r=a[i]/b[i-1]
        b[i]=b[i]-r*c[i-1]
        d[i]=d[i]-r*d[i-1]

    # /* The answers are stored in x[i] */
    d[n]=d[n]/b[n]

    for i in range (n-1 , 0 , -1):
        d[i]=(d[i]-c[i]*d[i+1]) / b[i]
    return

#++++++++++++++++++++++++++++++++++++++++
# n=99 aa=1.0 bb=2.0 alfa=0.5 bata=0.333333
#++++++++++++++++++++++++++++++++++++++++
n=99
aa=1.0
bb=2.0
alfa=0.5
bata=0.333333

a= [0.0 for i in range(0,n+1)]     #a [n] 矩陣
b= [0.0 for i in range(0,n+1)]     #a [n] 矩陣
c= [0.0 for i in range(0,n+1)]     #a [n] 矩陣
d= [0.0 for i in range(0,n+1)]     #a [n] 矩陣

y= [0.5 for i in range(0,n+1)]     #a [n] 矩陣


h=(bb-aa)/(n+1)

for k in range (1, 301):
    #x1=aa+h;
    x1=aa+h;
    b[1]=math.pow(h,2)* q ( (x1) , (y[1]) ) - 2.0
    c[1]=(1+(h/2.0)*p(  (x1) , (y[1]) ))
    d[1]=math.pow(h,2)*r( (x1) , (y[1]) ) - (1.0 - (h/2.0) * p(  (x1) , (y[1]))) * alfa

    for i in range (2 , n):
        x=aa+i*h;
        a[i]=1-(h/2.0)*p((x),(y[i]));
        b[i]=math.pow(h,2)*q((x),(y[i]))-2.0;
        c[i]=1+(h/2.0)*p((x),(y[i]));
        d[i]=math.pow(h,2)*r((x),(y[i]));

    xn=aa+n*h
    a[n]=1-(h/2.0)*p( (xn) , (y[n]) )
    b[n]=math.pow(h,2)*q((xn) , (y[n]))-2;
    d[n]=math.pow(h,2)*r((xn) , (y[n]))-(1+(h/2.0)*p( (xn) , (y[n]) ) )*bata;

    tridg(n,a,b,c,d)

    err=0.0
    for i in range (1 , n+1):
        err=err+abs(d[i]-y[i])
       
    if(err >0.001):
        for i in range (1, n+1):
            y[i]=d[i]
    else:
            break

#bound:
print("The iterations={%d}\n" %(k))
print(" x\t\ty(x)\t\t\tw(x)\t\t\t|y(x)-w(x)|\n")
print ("{%6.4f}\t\t{%10.7f}\t\t{%10.7f}\t\t{%10.7f}" %(aa,alfa, w(aa), abs(alfa-w(aa) ) ) )
for i in range (1 , n+1):
    x=aa+i*h;
    if(i%10==0):
        print ("{%6.4f}\t\t{%10.7f}\t\t{%10.7f}\t\t{%10.7f}" %(x,d[i],w(x),abs(d[i]-w(x) ) ))

print ("{%6.4f}\t\t{%10.7f}\t\t{%10.7f}\t\t{%10.7f}" %(bb,bata,w(bb),abs(bata-w(bb) ) ))

#for i in range (1,n+1) :      # //print the new matrix
 #   print( round(s1[i],4),"\t",end='')

輸出畫面
======== RESTART: F:\2018-09勤益科大數值分析\數值分析\PYTHON\EX7-4.py ===========
The iterations={5}

 x y(x) w(x) |y(x)-w(x)|

{1.0000} { 0.5000000} { 0.5000000} { 0.0000000}
{1.1000} { 0.4761904} { 0.4761905} { 0.0000000}
{1.2000} { 0.4545454} { 0.4545455} { 0.0000001}
{1.3000} { 0.4347825} { 0.4347826} { 0.0000001}
{1.4000} { 0.4166665} { 0.4166667} { 0.0000002}
{1.5000} { 0.3999998} { 0.4000000} { 0.0000002}
{1.6000} { 0.3846151} { 0.3846154} { 0.0000002}
{1.7000} { 0.3703701} { 0.3703704} { 0.0000003}
{1.8000} { 0.3571426} { 0.3571429} { 0.0000003}
{1.9000} { 0.3448273} { 0.3448276} { 0.0000003}
{2.0000} { 0.3333330} { 0.3333333} { 0.0000003}
>>>

數值分析 Numerical Analysis

網站著作權屬於國立臺灣師範大學 資訊工程學系
源自於
http://www.csie.ntnu.edu.tw/~u91029/

Search

 

About

About

Algorithm

Algorithm
Algorithm Design
» Simulation
» Backtracking
» Dynamic Programming
» cont. » cont. » cont.
Algorithm Analysis

 

Data Structure

Data » cont.
Set
Order
Sequence
Position
Region

Graph Theory

Graph » cont.
Tree » cont.
Directed Acyclic Graph
Bipartite Graph
Chordal Graph
Planar Graph
Component » cont.
Drawing

 

Combinatorial Optimization

Path » cont. » cont.
Spanning Tree » cont.
Flow » cont.
Cut » cont.
Matching » cont.
Circuit
Domination
Coloring
Labeling

Computational Geometry

Point » cont.
Shape » cont.
Polygon
Convex Hull » cont.
Triangulation
Neighbor

 

Computational Topology

Curve
Origami

Numerical Analysis

Function » cont. » cont.
Root Finding » cont.
Optimization » cont.
Interpolation
Regression » cont.
Classification
Transformation
Representation » cont.
Function Network

 

Signal Processing

Signal » Wave
Sequence » Poset
State » cont. » cont.
Function
Functional Equation

Number Theory

Number
Function » Permutation
Base
Divisor
Prime
Residue » Lattice
Game

 

Analytic Number Theory

Number
Sequence
Polynomial
Series

Stringology

String
String Searching
» cont. » cont. » cont.
String Alignment
» cont. » cont. » cont.
Language

 

Coding Theory

Code
Compression
Correction
Encryption
Transmission

Multimedia

Data
Text » Language
Audio » Audition
Image » Vision
Model » Motion
Graphics » Animation
Camera » Glass
Robot » Perception

 

Interdisciplinarity

Computer
Physics
Biology
Medicine
Economics
Education
Computer Science

例題2-2 1.15 ⋅ 10^ − 5 x ^2 + 1.4 ⋅ 10^ − 5 x − 1.962 ⋅ 10 ^− 2 = 0

例題2-2
利用 網站  https://www.mathway.com/zh/Precalculus

$1.15\cdot {10}^{-5}{x}^2+1.4\cdot {10}^{-5}x-1.962\cdot {10}^{-2}=0$

简化每一项。

点击了解更多的步骤...

$0.0000115x^2+0.000014x-0.01962=0$

从$0.0000115x^2+0.000014x-0.01962$中分解出因子$0.0000005$。

点击了解更多的步骤...

$0.0000005(23x^2+28x-39240)=0$

两边都除以$0.0000005$。$0$ 被任何非零数除结果都为$0$。

$(23x^2+28x-39240)=0$

将$23x^2+28x-39240$ 设置为等于$0$ 并且求解$x$。

点击了解更多的步骤...

$x=\frac{-14\pm 2\sqrt{225679}}{23}$

结果可用精确值和十进制两种方式表示。

Exact Form:

$x=\frac{-14+2\sqrt{225679}}{23},\frac{-14-2\sqrt{225679}}{23}$

Decimal Form:

$x=40.70059466\dots ,-41.91798597\dots$