Processing math: 100%

PPPK_01

262 days ago by fresl

Programiranje postupaka proračuna konstrukcija

 

Ishodi učenja predmeta:

nakon uspješno svladanog predmeta moći ćete

  • razumjeti i objasniti strukturu programa za proračun konstrukcija,
  • razumjeti, objasniti i predvidjeti posljedice neizbježnih aproksimacija u numeričkom modeliranju konstrukcija i ograničene (konačne) točnosti numeričkih proračuna,
  • oblikovati i napisati jednostavniji računalni program,
  • izmijeniti, prilagoditi i dograditi složeniji računalni program za proračun konstrukcija, dostupan u izvornom kodu,
  • surađivati u timu ili s timom koji oblikuje i piše složeni računalni program za proračun konstrukcija.

 

SageMath (Sage) — otvoreni programski sustav za simboličku matematiku (Computer Algebra System, CAS)

 

Literatura:

 

 

 

Elementi programiranja

uz malo matematičke i numeričke analize

 

 

Zadatak

 

                          (*)

                       

 

Treba odrediti minimum funkcije   L : x  (xxd)2+y2d+(xxv)2+y2v.

       
(x_d, y_d, x_v, y_v)
(x_d, y_d, x_v, y_v)
       
(xd,yd,xv,yv)
                                
                            

                                
       
       
x  (xxd)2+y2d+(xxv)2+y2v
                                
                            

                                
       
       
x  (x+1)2+16+(x5)2+1
                                
                            

                                
       
       

Matematička nas analiza uči da funkcija L poprima minimum u točki x u kojoj je, formalno, L(x)=0 ili, intuitivno–geometrijski, u točki u kojoj je tangenta na graf funkcije horizontalna, tako da funkcija u njoj ne raste niti ne pada. Dakle, ako je f=L, za nalaženje minimuma treba riješiti nelinearnu jednadžbu f(x)=0.

       
x  xxd(xxd)2+y2d+xxv(xxv)2+y2v
                                
                            

                                
       
x  x+1(x+1)2+16+x5(x5)2+1
                                
                            

                                

ili

       
x  x+1(x+1)2+16+x5(x5)2+1
                                
                            

                                
       
       
       
       
       
       
[x == (5*sqrt(x^2 + 2*x + 17) - sqrt(x^2 - 10*x + 26))/(sqrt(x^2 + 2*x + 17) +
sqrt(x^2 - 10*x + 26))]
[x == (5*sqrt(x^2 + 2*x + 17) - sqrt(x^2 - 10*x + 26))/(sqrt(x^2 + 2*x + 17) + sqrt(x^2 - 10*x + 26))]
       
x=5x2+2x+17x210x+26x2+2x+17+x210x+26
                                
                            

                                
       
[x == (19/5)]
[x == (19/5)]
       
3.80000000000000
3.80000000000000

Točno je rješenje   x=195=3,8.

 

 

Newton–Raphsonov postupak rješavanja nelinearnih jednadžbi f(x)=0  (za neku drugu, „opću” funkciju f)    

 

 


  • linearna aproksimacija funkcije f u točki x(k):    

ˉf(k)(x)=ax+b

a=f(x(k))

ˉf(k)(x(k))=ax(k)+b=f(x(k))f(x(k))x(k)+b=f(x(k))b=f(x(k))f(x(k))x(k)

ˉf(k)(x)=f(x(k))x+[f(x(k))f(x(k))x(k)]

 

  • nultočka funkcije ˉf(k):    

f(x(k))x+f(x(k))f(x(k))x(k)=0x=x(k)f(x(k))f(x(k))

 
  • iteracija:   x(k+1)=x(k)f(x(k))f(x(k))

 

  • završetak iteracije:     |f(x(k))|<τf   ili   k>N,   gdje je τf odabrana točnost, a N odabrani najveći broj koraka (za slučaj divergencije)

 

 

Još dva primjera:

 

 

 

 

 

       
x  (xxd)2((xxd)2+y2d)32(xxv)2((xxv)2+y2v)32+1(xxd)2+y2d+1(xxv)2+y2v
                                
                            

                                
       
x  (x+1)2((x+1)2+16)32(x5)2((x5)2+1)32+1(x+1)2+16+1(x5)2+1
                                
                            

                                

 

Programska funkcija definira se na sljedeći način:

    def naziv (parametri) :

         tijelo funkcije

Definicija funkcije započinje ključnom riječju def iza koje slijedi naziv funkcije. Dobro odabrani naziv trebao bi biti ključem značenja i svrhe funkcije. Nakon naziva se između okruglih zagrada navode parametri funkcije. Parametri su varijable koje pri upotrebi funkcije prihvaćaju vrijednosti koje se prenose u funkciju — koje „ulaze” u nju.  (Par zagrada treba napisati i ako funkcija nema parametara.)  Prvi redak definicije funkcije, koji se naziva i zaglavljem funkcije, završava dvotočkom.

Tijelo funkcije je niz naredaba koje propisuju kako funkcija radi ono čemu je namijenjena. Tijelo je programska realizacija algoritma.

Algoritam je razgovijetan, potpun i nedvosmislen opis postupka rješavanja određenoga zadatka u obliku niza jasno definiranih koraka. Svaki algoritam mora imati sljedeća svojstva:

  1. konačnost — izvođenje algoritma, neovisno o konkretnom skupu ulaznih/početnih podataka, mora završiti u konačnom broju koraka;
  2. jedinstveni početak — svako izvođenje algoritma mora početi istim korakom;
  3. jedinstveni nastavak — iza svakoga pojedinog koraka mora slijediti jednoznačno određeni korak;
  4. rješenje — svako izvođenje algoritma mora dovesti do rješenja ili završiti s naznakom da za zadani skup početnih podataka zadatak nije rješiv tim algoritmom.

Vrijednost koju funkcija izračunava i „vraća” naredbom return je vrijednost funkcije.

 

Upotreba funkcije naziva se i pozivom funkcije

       naziv (argumenti)

Nakon naziva se između okruglih zagrada navode funkcijski argumenti — vrijednosti koje se prenose u funkciju pridruživanjem parametrima.  (Argumenti se katkad nazivaju aktualnim parametrima, dok se za parametre tada upotrebljava naziv formalni parametri.)

  

       
x  x+1(x+1)2+16+x5(x5)2+1
                                
                            

                                
       
       
0.832050294337844
0.832050294337844
       
x  (x+1)2((x+1)2+16)32(x5)2((x5)2+1)32+1(x+1)2+16+1(x5)2+1
                                
                            

                                
       
       
1.04266924586348
1.04266924586348

 

Newton–Raphsonov postupak: 

x(0)=na neki način (primjerice bacanjem kocaka) pretpostavljeni broj

x(k+1)=x(k)f(x(k))f(x(k))

dok su   |f(x(k))|>τf   i   k<N    (oba uvjeta moraju biti zadovoljena)

 

Prva verzija programske funkcije newton_raphson():

       
       
3.80000000004635
3.80000000004635
       
1.52002854747479e-11
1.52002854747479e-11
       
NaN + NaN*I
NaN + NaN*I
       
       
0    5.00000000000000    0.832050294337844    1.04266924586348
1    4.20199977352474    0.168995782368261    0.534158108945978
2    3.88562194542862    0.0294799921364655    0.361472597627328
3    3.80406668417239    0.00133673130512302    0.329444455848224
4    3.80000915252663    3.00170821554424e-6    0.327966641072566
5    3.80000000004635    1.52002854747479e-11
3.80000000004635
0    5.00000000000000    0.832050294337844    1.04266924586348
1    4.20199977352474    0.168995782368261    0.534158108945978
2    3.88562194542862    0.0294799921364655    0.361472597627328
3    3.80406668417239    0.00133673130512302    0.329444455848224
4    3.80000915252663    3.00170821554424e-6    0.327966641072566
5    3.80000000004635    1.52002854747479e-11
3.80000000004635
       
0    7.00000000000000    1.78885438199983    0.111803398874990
1    -8.99999999999999    -1.89188589083065    0.0227223399272208
2    74.2610504415619    1.99848639858029    0.0000403829985181937
3    -49414.0501674566    -1.99999999651880    1.40901172836477e-13
4    1.41943459326446e13    2.00000000000000    0.000000000000000
5    -infinity    NaN    NaN
6    NaN    NaN + NaN*I    NaN + NaN*I
7    NaN + NaN*I    NaN + NaN*I    NaN + NaN*I
8    NaN + NaN*I    NaN + NaN*I    NaN + NaN*I
9    NaN + NaN*I    NaN + NaN*I    NaN + NaN*I
9    NaN + NaN*I    NaN + NaN*I
NaN + NaN*I
0    7.00000000000000    1.78885438199983    0.111803398874990
1    -8.99999999999999    -1.89188589083065    0.0227223399272208
2    74.2610504415619    1.99848639858029    0.0000403829985181937
3    -49414.0501674566    -1.99999999651880    1.40901172836477e-13
4    1.41943459326446e13    2.00000000000000    0.000000000000000
5    -infinity    NaN    NaN
6    NaN    NaN + NaN*I    NaN + NaN*I
7    NaN + NaN*I    NaN + NaN*I    NaN + NaN*I
8    NaN + NaN*I    NaN + NaN*I    NaN + NaN*I
9    NaN + NaN*I    NaN + NaN*I    NaN + NaN*I
9    NaN + NaN*I    NaN + NaN*I
NaN + NaN*I
       
x  x+1(x+1)2+16+x5(x5)2+1
                                
                            

                                
       
       
(-1, 4)
(5, 1)
(-1, 4)
(5, 1)

... u rješavanju nelinearnoga zadatka korisno je poznavati ga i razumjeti

 

Druga verzija:

       
       
3.80000000004635
3.80000000004635
       
-49414.0501674566
-49414.0501674566
       
-1.99999999651880
-1.99999999651880
       
       
(3.80000915252663, 5)
(3.80000915252663, 5)
       
(+Infinity, 4)
(+Infinity, 4)
       
Warning: maximal number of steps reached: 3
(3.80406668417239, 3)
Warning: maximal number of steps reached: 3
(3.80406668417239, 3)
       
3.80000000004635
1.52002854747479e-11
3.80000000004635
1.52002854747479e-11

 

Konačna verzija:

       
       
3.80000915252663
3.80000915252663
       
3.00170821554424e-6
3.00170821554424e-6
       
(3.80000915252663, 5)
(3.80000915252663, 5)
       
Traceback (click to the left of this block for traceback)
...
TypeError: unsupported operand parent(s) for +: '<type 'tuple'>' and
'Integer Ring'
Traceback (most recent call last):
  File "<stdin>", line 1, in <module>
  File "_sage_input_59.py", line 10, in <module>
    exec compile(u'open("___code___.py","w").write("# -*- coding: utf-8 -*-\\n" + _support_.preparse_worksheet_cell(base64.b64decode("ZnAgKHh4KQ=="),globals())+"\\n"); execfile(os.path.abspath("___code___.py"))
  File "", line 1, in <module>
    
  File "/tmp/tmpe1ew72/___code___.py", line 2, in <module>
    exec compile(u'fp (xx)
  File "", line 1, in <module>
    
  File "/tmp/tmpRn9J7H/___code___.py", line 4, in fp
    return (x + _sage_const_1 )/sqrt((x + _sage_const_1 )**_sage_const_2  + _sage_const_16 ) + (x - _sage_const_5 )/sqrt((x - _sage_const_5 )**_sage_const_2  + _sage_const_1 )
  File "sage/rings/integer.pyx", line 1792, in sage.rings.integer.Integer.__add__ (build/cythonized/sage/rings/integer.c:12197)
  File "sage/structure/coerce.pyx", line 1207, in sage.structure.coerce.CoercionModel.bin_op (build/cythonized/sage/structure/coerce.c:10896)
TypeError: unsupported operand parent(s) for +: '<type 'tuple'>' and 'Integer Ring'
       
3.00170821554424e-6
3.00170821554424e-6
       
3.80000000004635
6
3.80000000004635
6
       
1.52002854747479e-11
1.52002854747479e-11
       
Traceback (click to the left of this block for traceback)
...
TypeError: 'sage.rings.real_mpfr.RealNumber' object is not iterable
Traceback (most recent call last):
  File "<stdin>", line 1, in <module>
  File "_sage_input_63.py", line 10, in <module>
    exec compile(u'open("___code___.py","w").write("# -*- coding: utf-8 -*-\\n" + _support_.preparse_worksheet_cell(base64.b64decode("eHgsIHN0ID0gbmV3dG9uX3JhcGhzb24gKGZwLCBkZnAsIDUuLCB0b2xfZiA9IDEuZS03KQ=="),globals())+"\\n"); execfile(os.path.abspath("___code___.py"))
  File "", line 1, in <module>
    
  File "/tmp/tmpYk2ycc/___code___.py", line 3, in <module>
    exec compile(u'xx, st = newton_raphson (fp, dfp, _sage_const_5p , tol_f = _sage_const_1pen7 )
  File "", line 1, in <module>
    
TypeError: 'sage.rings.real_mpfr.RealNumber' object is not iterable
       
3.80000000004635
3.80000000004635
       
step = 0   x = 5.00000000000000   f(x) = 0.832050294337844
step = 1   x = 4.20199977352474   f(x) = 0.168995782368261
step = 2   x = 3.88562194542862   f(x) = 0.0294799921364655
step = 3   x = 3.80406668417239   f(x) = 0.00133673130512302
step = 4   x = 3.80000915252663   f(x) = 3.00170821554424e-6
step = 5   x = 3.80000000004635   f(x) = 1.52002854747479e-11
step = 6   x = 3.80000000000000   f(x) = 0.000000000000000
step = 0   x = 5.00000000000000   f(x) = 0.832050294337844
step = 1   x = 4.20199977352474   f(x) = 0.168995782368261
step = 2   x = 3.88562194542862   f(x) = 0.0294799921364655
step = 3   x = 3.80406668417239   f(x) = 0.00133673130512302
step = 4   x = 3.80000915252663   f(x) = 3.00170821554424e-6
step = 5   x = 3.80000000004635   f(x) = 1.52002854747479e-11
step = 6   x = 3.80000000000000   f(x) = 0.000000000000000
       
step = 0   x = 5.00000000000000   f(x) = 0.832050294337844
step = 1   x = 4.20199977352474   f(x) = 0.168995782368261
step = 2   x = 3.88562194542862   f(x) = 0.0294799921364655
step = 3   x = 3.80406668417239   f(x) = 0.00133673130512302
Warning: maximal number of steps reached: 4
3.80000915252663
step = 0   x = 5.00000000000000   f(x) = 0.832050294337844
step = 1   x = 4.20199977352474   f(x) = 0.168995782368261
step = 2   x = 3.88562194542862   f(x) = 0.0294799921364655
step = 3   x = 3.80406668417239   f(x) = 0.00133673130512302
Warning: maximal number of steps reached: 4
3.80000915252663
       
step = 0   x = 3.78900000000000   f(x) = -0.00358578168341916
step = 1   x = 3.80006689346948   f(x) = 0.0000219394162994657
step = 2   x = 3.80000000247579   f(x) = 8.11967382219336e-10
3.80000000247579
8.11967382219336e-10
step = 0   x = 3.78900000000000   f(x) = -0.00358578168341916
step = 1   x = 3.80006689346948   f(x) = 0.0000219394162994657
step = 2   x = 3.80000000247579   f(x) = 8.11967382219336e-10
3.80000000247579
8.11967382219336e-10
       
3.80000000000000
0.000000000000000
7
3.80000000000000
0.000000000000000
7
       
3.80000000010167
3.80000000010167
       
Tangent slope f'(x) = 1.04436621297613e-16  <  1.00000000000000e-7
+Infinity
Tangent slope f'(x) = 1.04436621297613e-16  <  1.00000000000000e-7
+Infinity
       
Tangent slope f'(x) = 1.73422661322621e-12  <  1.00000000000000e-7
+Infinity
Tangent slope f'(x) = 1.73422661322621e-12  <  1.00000000000000e-7
+Infinity
       
3.80000006304201
3.80000006304201
       
Tangent slope f'(x) = 3.77016804277777e-14  <  1.00000000000000e-7
+Infinity
Tangent slope f'(x) = 3.77016804277777e-14  <  1.00000000000000e-7
+Infinity

 

Jedno poopćenje:

       
x  xxd(xxd)2+y2d+xxv(xxv)2+y2v
                                
                            

                                
       
x  x+1(x+1)2+16+x5(x5)2+1
                                
                            

                                
       
       
(-1, 4)
(5, 1)
(-1, 4)
(5, 1)
       
0.832050294337844
0.832050294337844
0.832050294337844
0.832050294337844
       
x  (xxd)2((xxd)2+y2d)32(xxv)2((xxv)2+y2v)32+1(xxd)2+y2d+1(xxv)2+y2v
                                
                            

                                
       
       
Traceback (click to the left of this block for traceback)
...
TypeError: gp() takes exactly 3 arguments (1 given)
Traceback (most recent call last):
  File "<stdin>", line 1, in <module>
  File "_sage_input_81.py", line 10, in <module>
    exec compile(u'open("___code___.py","w").write("# -*- coding: utf-8 -*-\\n" + _support_.preparse_worksheet_cell(base64.b64decode("bmV3dG9uX3JhcGhzb24gKGdwLCBkZ3AsIDUuKQ=="),globals())+"\\n"); execfile(os.path.abspath("___code___.py"))
  File "", line 1, in <module>
    
  File "/tmp/tmpfrMX8w/___code___.py", line 3, in <module>
    exec compile(u'newton_raphson (gp, dgp, _sage_const_5p )
  File "", line 1, in <module>
    
  File "/tmp/tmpQlTEPJ/___code___.py", line 9, in newton_raphson
    fk = f (xk)
TypeError: gp() takes exactly 3 arguments (1 given)
       
Traceback (click to the left of this block for traceback)
...
TypeError: gp() takes exactly 3 arguments (1 given)
Traceback (most recent call last):
  File "<stdin>", line 1, in <module>
  File "_sage_input_82.py", line 10, in <module>
    exec compile(u'open("___code___.py","w").write("# -*- coding: utf-8 -*-\\n" + _support_.preparse_worksheet_cell(base64.b64decode("bmV3dG9uX3JhcGhzb24gKGdwLCBkZ3AsIDUuLCB4eWQsIHh5dik="),globals())+"\\n"); execfile(os.path.abspath("___code___.py"))
  File "", line 1, in <module>
    
  File "/tmp/tmpvnTklM/___code___.py", line 3, in <module>
    exec compile(u'newton_raphson (gp, dgp, _sage_const_5p , xyd, xyv)
  File "", line 1, in <module>
    
  File "/tmp/tmpQlTEPJ/___code___.py", line 9, in newton_raphson
    fk = f (xk)
TypeError: gp() takes exactly 3 arguments (1 given)
       
       
0.832050294337844
0.832050294337844
       
<function f_ at 0x7f562dbd0c08>
<function f_ at 0x7f562dbd0c08>
       
<function f_ at 0x7f562dbd0b18>
<function f_ at 0x7f562dbd0b18>
       
0.832050294337844
0.832050294337844
       
       
<function <lambda> at 0x7f562de13668>
<function <lambda> at 0x7f562de13668>
       
0.832050294337844
0.832050294337844
       
3.80000915252663
3.80000915252663

 

Minimizacija primjenom Newton-Raphsonova postupka:

       
       
step = 0   x = 6.00000000000000   f(x) = 1.57534992331101
step = 1   x = 1.89843432768716   f(x) = -0.364994091915365
step = 2   x = 4.15652661816414   f(x) = 0.145393647790519
step = 3   x = 3.86816299760975   f(x) = 0.0232325785022762
step = 4   x = 3.80257728615178   f(x) = 0.000846462445067031
step = 5   x = 3.80000367572973   f(x) = 1.20550697590982e-6
step = 6   x = 3.80000000000748   f(x) = 2.45159448297727e-12
3.80000000000748
7.81024967590665
step = 0   x = 6.00000000000000   f(x) = 1.57534992331101
step = 1   x = 1.89843432768716   f(x) = -0.364994091915365
step = 2   x = 4.15652661816414   f(x) = 0.145393647790519
step = 3   x = 3.86816299760975   f(x) = 0.0232325785022762
step = 4   x = 3.80257728615178   f(x) = 0.000846462445067031
step = 5   x = 3.80000367572973   f(x) = 1.20550697590982e-6
step = 6   x = 3.80000000000748   f(x) = 2.45159448297727e-12
3.80000000000748
7.81024967590665
       
       
False
True
False
True
       
False
False
True
False
False
True
       
       
(3.80000000000000, 7)
(3.80000000000000, 7)
       
False
True
False
True
       
x |--> sqrt((x - x_d)^2 + y_d^2) + sqrt((x - x_v)^2 + y_v^2)
x |--> sqrt((x - x_d)^2 + y_d^2) + sqrt((x - x_v)^2 + y_v^2)
       
3.80000000000000
3.80000000000000


(*)  Crteži Dona Quijotea, Sancha Panze i vjetrenjača posuđeni su s omotnice albuma Windmill Tilter Kena Wheelera i Orkestra Johna Dankwortha iz 1969. godine.