Диссертация (1144294), страница 20
Текст из файла (страница 20)
mat ’] , ’ currPointsPS0 ’);preIterPoincareMap = 1;elseif isCont == 1lastIndex = getLastFigIndex ( DIR ) ;load ([ DIR ’ currPointsPS0 . mat ’] , ’ currPointsPS0 ’);preIterPoincareMap = 0;elsereturn ;end151152153154% % Main loopfor iPoincareMap = 1 : numPoincareMaps - preIterPoincareMap142155[ currPointsBackPS0 , currPointsPS1_1 , currPointsPS1_2 ] =execPoincareMap ( currPointsPS0 , odeLorenzLike , [0 T ],optionsOnPS0 , optionsOnPS1 ) ;156currIndex = lastIndex + iPoincareMap ;fileName = [ DIR num2str ( currIndex , ’ %02 d ’) ];157158159160plotPoincareSec ( equilPoints , vS , vSS , rect3D , sepa , frVertex ,currPointsPS0 , currPointsPS1_1 , currPointsPS1_2 ,currPointsBackPS0 , axisSize3D , fileName );161currPointsPS0 = currPointsBackPS0 ;save ([ DIR ’ currPointsPS0 .
mat ’] , ’ currPointsPS0 ’);162163164end165166warning ( ’ on ’ , ’ MATLAB : odearguments : RelTolIncrease ’) ;167168endЛистинг В.17:execPoincareMap.m – функция, выполняющая отображениеПуанкаре.12345function [ pointsBackOnPS0 , pointsOnPS1_1 , pointsOnPS1_2 ] =execPoincareMap ( pointsOnPS0 , ode , tInt , optionsOnPS0 ,optionsOnPS1 )if ~ isempty ( pointsOnPS0 )numWorkers = feature ( ’ numCores ’) ;poolObj = parpool ( numWorkers ) ;cleaner = onCleanup (@ () delete ( poolObj ) ) ;6789101112131415[ rPS0 , cPS0 ] = size ( pointsOnPS0 ) ;startPointPS0 = pointsOnPS0 (: , 1:3) ;colorPS0 = pointsOnPS0 (: , 4:6) ;pointsOnPS1 = zeros ( rPS0 , cPS0 +1) ;parfor iPS0 = 1 : rPS0warning ( ’ off ’ , ’ MATLAB : odearguments : RelTolIncrease ’) ;[ ~ , currTraj , ~ , ~ , currIe ] = ode45 ( ode , tInt ,startPointPS0 ( iPS0 , :) , optionsOnPS1 ) ;pointsOnPS1 ( iPS0 , :) = [ currTraj ( end , :) , colorPS0 ( iPS0, :) , currIe ];warning ( ’ on ’ , ’ MATLAB : odearguments : RelTolIncrease ’) ;143161718endpointsOnPS1_1 = pointsOnPS1 ( pointsOnPS1 (: , end ) == 1 , 1 :end -1) ;pointsOnPS1_2 = pointsOnPS1 ( pointsOnPS1 (: , end ) == 2 , 1 :end -1) ;19if ~ isempty ( pointsOnPS1_1 )[ rPS1_1 , cPS1_1 ] = size ( pointsOnPS1_1 ) ;pointsBackOnPS0_1 = zeros ( rPS1_1 , cPS1_1 ) ;startPointPS1_1 = pointsOnPS1_1 (: , 1:3) ;colorPS1_1 = pointsOnPS1_1 (: , 4:6) ;parfor iPS1_1 = 1 : rPS1_1warning ( ’ off ’ , ’ MATLAB : odearguments : RelTolIncrease ’20212223242526);2728293031323334[ ~ , currTraj ] = ode45 ( ode , tInt , startPointPS1_1 (iPS1_1 , :) , optionsOnPS0 ) ;pointsBackOnPS0_1 ( iPS1_1 , :) = [ currTraj ( end , :) ,colorPS1_1 ( iPS1_1 , :) ];warning ( ’ on ’, ’ MATLAB : odearguments : RelTolIncrease ’) ;endelsepointsBackOnPS0_1 = [];endpointsBackOnPS0 = pointsBackOnPS0_1 ;35if ~ isempty ( pointsOnPS1_2 )[ rPS1_2 , cPS1_2 ] = size ( pointsOnPS1_2 ) ;pointsBackOnPS0_2 = zeros ( rPS1_2 , cPS1_2 ) ;startPointPS1_2 = pointsOnPS1_2 (: , 1:3) ;colorPS1_2 = pointsOnPS1_2 (: , 4:6) ;parfor iPS1_2 = 1 : rPS1_2warning ( ’ off ’ , ’ MATLAB : odearguments : RelTolIncrease ’36373839404142);4344454647[ ~ , currTraj ] = ode45 ( ode , tInt , startPointPS1_2 (iPS1_2 , :) , optionsOnPS0 ) ;pointsBackOnPS0_2 ( iPS1_2 , :) = [ currTraj ( end , :) ,colorPS1_2 ( iPS1_2 , :) ];warning ( ’ on ’, ’ MATLAB : odearguments : RelTolIncrease ’) ;endelse144pointsBackOnPS0_2 = [];48endpointsBackOnPS0 = [ pointsBackOnPS0 ; pointsBackOnPS0_2 ];4950else515253545556errorMessage = sprintf ( ’ Error : No points on Poincaresection \ n ’) ;uiwait ( warndlg ( errorMessage ) );return ;endendЛистинг В.18:rectangle3d.m – функция, задающая вершины прямоугольнойрамки на сечении Пуанкаре Σin .1function rect3d = rectangle3d (P0 , vHor , vVert , rectHorSize ,rectVertSize )2p1p2p3p43456====P0P0P0P0++rectHorSizerectHorSizerectHorSizerectHorSize****vHorvHorvHorvHor++rectVertSizerectVertSizerectVertSizerectVertSize****vVert ;vVert ;vVert ;vVert ;7rect3d = [ p1 (:) , p2 (:) , p3 (:) , p4 (:) ];89endЛистинг В.19:frameVertex.m – функция, задающая вершины полурамки насечении Пуанкаре Σin .12function frVertex = frameVertex ( midPoint , vHor , vVert , ...frThickness , frInHalfWidth , frHalfHeight )34567891011P2P1P8P3P4P5P6P7========midPoint + frInHalfWidth * vHor ;P2 + frThickness * vHor ;P1 + frHalfHeight * vVert ;P2 + ( frHalfHeight - frThickness ) * vVert ;P3 - 2 * frInHalfWidth * vHor ;midPoint - frInHalfWidth * vHor ;P5 - frThickness * vHor ;P6 + frHalfHeight * vVert ;1213frVertex = [ P1 ; P2 ; P3 ; P4 ; P5 ; P6 ; P7 ; P8 ];1451415endЛистинг В.20:frameGrid.m – функция, генерирующая начальную сеткуточек в виде полурамки на сечении Пуанкаре Σin .12function grid = frameGrid ( midPoint , vHor , vVert , ...frThickness , frInHalfWidth , frHalfHeight , numPoints )3456grid = []; currIndex = 1;% % 1 st rectanglepointRect1 = midPoint + frInHalfWidth * vHor ;78910stepSize = frThickness / ( numPoints -1) ;numPointsVert = floor ( frHalfHeight / stepSize ) + 1;col = jet ( numPointsVert ) ;1112131415161718for iPointX = 1 : numPointsxp = pointRect1 + ( iPointX -1) * stepSize * vHor ;for iPointY = 1 : numPointsVertgrid = [ grid ; xp + ( iPointY -1) * stepSize * vVert , col (iPointY , :) ];currIndex = currIndex +1;endend192021% % 2 nd rectanglepointRect2 = midPoint - frInHalfWidth * vHor ;2223242526272829for iPointX = 1 : numPointsxp = pointRect2 - ( iPointX -1) * stepSize * vHor ;for iPointY = 1 : numPointsVertgrid = [ grid ; xp + ( iPointY -1) * stepSize * vVert , col (iPointY , :) ];currIndex = currIndex +1;endend3031numPointsHor = 2 * frInHalfWidth / stepSize + 1;3233% % 3 nd rectangle14634pointRect3 = midPoint + frInHalfWidth * vHor + ( numPointsVert -1) *stepSize * vVert ;3536373839404142for iPointX = 1 : numPointsHor -1xp = pointRect3 - iPointX * stepSize * vHor ;for iPointY = 1 : numPointsgrid = [ grid ; xp - ( iPointY -1) * stepSize * vVert , col ( end- ( iPointY -1) ,:) ];currIndex = currIndex +1;endend434445grid (: , 1:3) = grid (: , 1:3) + stepSize * vVert ;endЛистинг В.21:plotPoincareSec.m – функция для отрисовки сетки точек насечении Пуанкаре при последовательном применении отображения Пуанкареи сохранения соответствующих картинок.1function plotPoincareSec ( equilibria , vS , vSS , rect , separatrix ,frVertex , pointsPS0 , pointsPS1_1 , pointsPS1_2 , pointsBackPS0 ,axisSize3D , fileName )234S0 = equilibria (1 , :) ; S1 = equilibria (2 , :) ; S2 = equilibria (3 , :);xPS0 = rect (1 , :) ; yPS0 = rect (2 , :) ; zPS0 = rect (3 , :) ;5678fx = @ (s , t ) vS (1) * s + vSS (1) * t ;fv = @ (s , t ) vS (2) * s + vSS (2) * t ;fu = @ (s , t ) vS (3) * s + vSS (3) * t ;910epsManifold = 1 e1 ;111213hFig1 = figure ( ’ Visible ’, ’ off ’) ; hold on ;axFig1 = gca ;141516% % Equilibria :plot3 ( axFig1 , S0 (1) , S0 (2) , S0 (3) , ’.
’ , ’ markersize ’ , 20 , ’ Color ’ ,’ black ’) ;1471718plot3 ( axFig1 , S1 (1) , S1 (2) , S1 (3) , ’. ’ , ’ markersize ’ , 20 , ’ Color ’ ,’ black ’) ;plot3 ( axFig1 , S2 (1) , S2 (2) , S2 (3) , ’. ’ , ’ markersize ’ , 20 , ’ Color ’ ,’ black ’) ;192021% % 2 D stable manifoldfsurf ( axFig1 , fx , fv , fu , [ - epsManifold , epsManifold ] , ’ EdgeColor ’ ,’ none ’, ’ FaceColor ’ , [0.9 , 0.9 , 0.9]) ;222324252627% % Separatrixif ~ isempty ( separatrix )plot3 ( axFig1 , separatrix (: ,1) , separatrix (: ,2) , separatrix (: ,3));plot3 ( axFig1 , - separatrix (: ,1) , - separatrix (: ,2) , separatrix(: ,3) ) ;end282930313233343536% % Poincare section :hPS1 = fill3 ( axFig1 , xPS0 , yPS0 , zPS0 , ’r ’) ; set ( hPS1 , ’ FaceAlpha ’ ,0.1) ;plotDashedLine ( axFig1 , [( xPS0 (1) + xPS0 (4) ) / 2 , ...( yPS0 (1) + yPS0 (4) ) / 2 , ...( zPS0 (1) + zPS0 (4) ) / 2] , ...[( xPS0 (2) + xPS0 (3) ) / 2 , ...( yPS0 (2) + yPS0 (3) ) / 2 , ...( zPS0 (2) + zPS0 (3) ) / 2]) ;37383940414243444546% % Half - frameplotDashedLine ( axFig1 ,plotDashedLine ( axFig1 ,plotDashedLine ( axFig1 ,plotDashedLine ( axFig1 ,plotDashedLine ( axFig1 ,plotDashedLine ( axFig1 ,plotDashedLine ( axFig1 ,plotDashedLine ( axFig1 ,frVertex (1 ,frVertex (2 ,frVertex (3 ,frVertex (4 ,frVertex (5 ,frVertex (6 ,frVertex (7 ,frVertex (8 ,:) ,:) ,:) ,:) ,:) ,:) ,:) ,:) ,frVertex (2 ,frVertex (3 ,frVertex (4 ,frVertex (5 ,frVertex (6 ,frVertex (7 ,frVertex (8 ,frVertex (1 ,:) ) ;:) ) ;:) ) ;:) ) ;:) ) ;:) ) ;:) ) ;:) ) ;474849% % Grid of pointsscatter3 ( axFig1 , pointsPS0 (: ,1) , pointsPS0 (: ,2) , pointsPS0 (: ,3) , 36 ,pointsPS0 (: ,4:6) , ’ filled ’) ;14850515253545556if ~ isempty ( pointsPS1_1 )scatter3 ( axFig1 , pointsPS1_1 (: ,1) , pointsPS1_1 (: ,2) , pointsPS1_1(: ,3) , 36 , pointsPS1_1 (: ,4:6) , ’ filled ’) ;endif ~ isempty ( pointsPS1_2 )scatter3 ( axFig1 , pointsPS1_2 (: ,1) , pointsPS1_2 (: ,2) , pointsPS1_2(: ,3) , 36 , pointsPS1_2 (: ,4:6) , ’ filled ’) ;end57585960616263grid ( axFig1 , ’on ’) ;if ~ isempty ( axisSize3D )axis ( axFig1 , axisSize3D ) ;elseaxis ( axFig1 , ’ auto ’) ;end6465666768xlabel ( axFig1 , ’x ’) ; ylabel ( axFig1 , ’y ’) ; zlabel ( axFig1 , ’z ’) ;% view ( ax1 , -104 ,14) ;view ( axFig1 , 3) ;hold off ;6970717273set ( hFig1 , ’ Position ’ , [1 , 1 , 500 , 500]) ;% print ( hFig2 , ’- r300 ’ , [ fileName , ’.
eps ’] , ’- depsc ’ , ’- tiff ’) ;savefig ( hFig1 , [ fileName , ’_ ( PS0 - > PS1 ) . fig ’] , ’ compact ’) ;close ( hFig1 ) ;74757677%%hFig2 = figure ( ’ Visible ’, ’ off ’) ; hold on ;axFig2 = gca ;78798081plot3 ( axFig2 , S0 (1) , S0 (2) , S0 (3) , ’. ’ , ’ markersize ’ , 20 , ’ Color ’ ,’ black ’) ;plot3 ( axFig2 , S1 (1) , S1 (2) , S1 (3) , ’. ’ , ’ markersize ’ , 20 , ’ Color ’ ,’ black ’) ;plot3 ( axFig2 , S2 (1) , S2 (2) , S2 (3) , ’. ’ , ’ markersize ’ , 20 , ’ Color ’ ,’ black ’) ;8283% % 2 D stable manifold14984fsurf ( axFig2 , fx , fv , fu , [ - epsManifold , epsManifold ] , ’ EdgeColor ’ ,’ none ’, ’ FaceColor ’ , [0.9 , 0.9 , 0.9]) ;8586878889if ~ isempty ( separatrix )plot3 ( axFig2 , separatrix (: ,1) , separatrix (: ,2) , separatrix (: ,3));plot3 ( axFig2 , - separatrix (: ,1) , - separatrix (: ,2) , separatrix(: ,3) ) ;end9091929394959697hPS1 = fill3 ( axFig2 , xPS0 , yPS0 , zPS0 , ’r ’) ; set ( hPS1 , ’ FaceAlpha ’ ,0.1) ;plotDashedLine ( axFig2 , [( xPS0 (1) + xPS0 (4) ) / 2 , ...( yPS0 (1) + yPS0 (4) ) / 2 , ...( zPS0 (1) + zPS0 (4) ) / 2] , ...[( xPS0 (2) + xPS0 (3) ) / 2 , ...( yPS0 (2) + yPS0 (3) ) / 2 , ...( zPS0 (2) + zPS0 (3) ) / 2]) ;9899100101102103104105106107% % Half - frameplotDashedLine ( axFig2 ,plotDashedLine ( axFig2 ,plotDashedLine ( axFig2 ,plotDashedLine ( axFig2 ,plotDashedLine ( axFig2 ,plotDashedLine ( axFig2 ,plotDashedLine ( axFig2 ,plotDashedLine ( axFig2 ,frVertex (1 ,frVertex (2 ,frVertex (3 ,frVertex (4 ,frVertex (5 ,frVertex (6 ,frVertex (7 ,frVertex (8 ,:) ,:) ,:) ,:) ,:) ,:) ,:) ,:) ,frVertex (2 ,frVertex (3 ,frVertex (4 ,frVertex (5 ,frVertex (6 ,frVertex (7 ,frVertex (8 ,frVertex (1 ,:) ) ;:) ) ;:) ) ;:) ) ;:) ) ;:) ) ;:) ) ;:) ) ;108109110111112113114115116if ~ isempty ( pointsPS1_1 )scatter3 ( axFig2 , pointsPS1_1 (: ,1) , pointsPS1_1 (: ,2) , pointsPS1_1(: ,3) , 36 , pointsPS1_1 (: ,4:6) , ’ filled ’) ;endif ~ isempty ( pointsPS1_2 )scatter3 ( axFig2 , pointsPS1_2 (: ,1) , pointsPS1_2 (: ,2) , pointsPS1_2(: ,3) , 36 , pointsPS1_2 (: ,4:6) , ’ filled ’) ;endscatter3 ( axFig2 , pointsBackPS0 (: ,1) , pointsBackPS0 (: ,2) ,pointsBackPS0 (: ,3) , 36 , pointsBackPS0 (: ,4:6) , ’ filled ’) ;grid ( axFig2 , ’on ’) ;150117118119120121122123124if ~ isempty ( axisSize3D )axis ( axFig2 , axisSize3D ) ;elseaxis ( axFig2 , ’ auto ’) ;endxlabel ( axFig2 , ’x ’) ; ylabel ( axFig2 , ’v ’) ; zlabel ( axFig2 , ’u ’) ;view ( axFig2 , 3) ;hold off ;125126127128129130set ( hFig2 , ’ Position ’ , [1 , 1 , 500 , 500]) ;% print ( hFig2 , ’- r300 ’ , [ fileName , ’.
eps ’] , ’- depsc ’ , ’- tiff ’) ;savefig ( hFig2 , [ fileName , ’_ ( PS1 - > PS0 ) . fig ’] , ’ compact ’) ;close ( hFig2 ) ;endЛистинг В.22:plotDashedLine.m – вспомогательная функция для отрисовкипунктирной линии.1function h = plotDashedLine ( ax , pBegin , pEnd )23h = plot3 ( ax , [ pBegin (1) , pEnd (1) ] , [ pBegin (2) , pEnd (2) ] , [ pBegin(3) , pEnd (3) ], ’ -- ’ , ’ LineWidth ’ , 1.5 , ’ Color ’ , ’ black ’) ;45endЛистинг В.23:getLastFigIndex.m – вспомогательная функция дляопределения номера последней итерации отображения Пуанкаре исоответствующей картинки сечения Пуанкаре в текущей папке.1function lastIndex = getLastFigIndex ( DIR )2figFiles = dir ([ DIR ’/* _ ( PS0 - > PS1 ) .
fig ’ ]) ;34if isempty ( figFiles )lastIndex = -1;elselastFigName = figFiles ( end ) . name ;lastFigNameSplited = strsplit ( lastFigName , ’_ ’) ;lastIndex = str2double ( cell2mat ( lastFigNameSplited (1) ) ) ;end56789101112end.
![](https://s.studizba.com/z.php?f=/templates/si/i/noavatar.png&w=100&h=100&t=1"){kind=avatar}
