problème pour traduire un code Matlab en Mathematica
Répondre à la discussion
Affichage des résultats 1 à 4 sur 4

problème pour traduire un code Matlab en Mathematica



  1. #1
    invitea41a00b4

    problème pour traduire un code Matlab en Mathematica


    ------

    Bonjour,
    J'ai recu plusieurs codes écrits en Matlab, et je voudrais les convertir en Mathematica. Je ne sais pas me servir de Mathlab et mon établissement ne le possede pas donc je dois absolument les convertir en Mathematica . Est ce que quelqu'un peut m'aider ? Merci d'avance
    Voici les codes:

    (1er code)
    Code:
    (Defining the constant surface for the Poincare section
     Zerocross2.m)
     
    function [v,ist,df] = zerocross2(t,x,varargin)% varargin absorbs unwanted
     parameters flag,g,l1, etc
     % event function for test_poincare.m
     v = sin(x(2)); % event is x(2)=0 (ie when v is zero)
     %v = x(2);
     ist = 1; % if true, terminate evolution when this event occurs
     df = 1; % increasing sense only
     
    (2ème code)(Producing the trajectory of the outer bob, Poincare section, angular velocity graph and energy of the system
     PoincareTrajecAngleEnergy.m) 
    
    clear
     %Code mainly adapted from Alexander Barnett
     %declare time step
     T = 10;
     %y0 = [0;0;-2;0]; % close to regular
     %y0 = [pi;pi;0;0];
     % y0 = [pi;pi;.5;0]; % chaotic
     % y0=[0.5233;0;0.5233;0]%periodic
     y0 = [0.2,0.2828,0,0]%perfect periodic with energy = 0.7809
     %y0 = [0.2,-0.2828,0,0]%perfect quasiperiodic with energy = 0.7809
     %y0 = [27.8;0;1.22;2.62] %perfect KAM scenario
     %y0 = [0,0,2,20]
     %y0 = [0,0,1,20];
     y0 =[pi;0;.5;.5];%high energy
     s=ode45(@doublependulum2,[0,T],y0,[],0,1,1,1,1);
     t = 0:0.01:T; x = deval(s,t)';
     x1=sin(x(:,1));
     y1=-cos(x(:,1));
     x2=x1+sin(x(:,2));
     y2=y1-cos(x(:,2));
     figure; plot(x2,y2); axis equal; xlabel('position');ylabel('vel ocity');title
     ('trajectory')
     tmax = 1e2; % max time to wait until next intersection
     ns = 1000; % how many intersection (iterations of P map)
     tp = nan*(1:ns); yp = nan*zeros(numel(y0),ns); yi = y0; % init arrays
     figure;
     for n=1:ns
     s = ode45(@doublependulum2, [0, tmax], yi,
     odeset('Events',@zerocross2,'a bstol',1e-9),0,9.8,1,1,1,1);
     if isempty(s.xe), disp('no intersection found!');
     else tp(n) = s.xe(end); yi = s.ye(:,end); yp(:,n) = yi;
     %plot3(mod(yi(1),2*pi),yi(3),y i(4),'+'); hold on; if mod(n,10)==0, drawnow; end
     % note it's in 3d now
     plot(yi(1),yi(4));hold on; if mod(n,10)==0, drawnow; end
     end
     end
     xlabel('x_1'); ylabel('x_3'); zlabel('x_4'); axis vis3d;
     title('Poincare section (sin(x_2)=0, increasing)');
     figure;
     plot(t,x(:,3),'magenta',t,x(:, 4),'black');title('Angular velocities');
     energy = hamiltonian(y0)
     
    (3ème code)(Calculates the lyapunov exponent, the average lyapunov exponent and traces the separation of close trajectories
     Crudelyap.m)
     
    clear
     nruns =100;
     %for calculating average later
     h = zeros(1,nruns);
     %indexing of h
     k = 1;
     for j=1:nruns
     % y0 = [pi+0.1*rand(1);pi;.5;0];
     % y0 = [0+0.1*rand(1);0;-2;0];
     %y0 = [pi;pi;10.515;-2.17];%periodic with T = 0.5
     %y0 = [];%quasiperiodic;
     y0 = [pi+0.1*rand(1);pi;0.5;0];%chaotic
     %y0 = [pi;pi;0;0];
     %y0 = [0;0;-2;0];%close to regular
     %y0 = [27.8;0;1.22;2.62];
     %y0 = [pi;0;0.5;0.5];
     %y0 = [0;0;2;20];
     %y0 = [pi;0;5;5];%high energy
     %y0 = [0.2;0.2828;0;0];%perfect periodic with energy = 0.7809
     %y0 = [0.2+0.1*rand(1);-0.2828;0;0];
     %y0 = [0.2;-0.2828;0;0];%perfect quasiperiodic with energy = 0.7809
     T = 20;
     o = odeset('abstol',1e-14);
     s=ode45(@doublependulum2,[0,T],y0,o,0,1,1,1,1); % note changes
     eps = 1e-9;
     s2 = ode45(@doublependulum2,[0,T],y0+[eps;0;0;0],o,0,1,1,1,1); % note changes
     t = 0:0.01:T; x = deval(s,t)';xe = deval(s2,t)';
     x1=sin(x(:,1));
     y1=-cos(x(:,1));
     x2=x1+sin(x(:,2));
     y2=y1-cos(x(:,2));
     xe1 = sin(xe(:,1));
     ye1=-cos(xe(:,1));
     xe2 = xe1+sin(xe(:,2));
     ye2= ye1-cos(xe(:,2));
     %store the size
     N = numel(ye2);
     %figure; plot(t, [y2 ye2], '-');title('Trajectories of two intial points with a
     distance 10^-9')
     logsep = zeros(1,N);
     %loop fills matrix logsep with natural log of the changing separation
     %between the initial trajectories
     for i = 1:N
     logsep(i) = (log(abs(y2(i)-ye2(i)))-log(eps));
     end
     %figure; plot(t,logsep);title('lyapunov exponent');
     t1=5;
     t2=0.7;
     %stores the exponent in a matric
     h(k) = (logsep(t1/0.01)-logsep(t2/0.01))/(t1-t2);
     k=k+1;
     end
     figure;
     %average lyapunov exponent
     plot(1:j,h)
     title('average lyapunov exponent');
     
    (4ème code)
     
    clear
     T = 10;
     %y0 = [0;0;0;10]; %really good
     %y0 = [pi;pi;0.5;0];
     %y0 = [0;0;-2;0];
     %y0=[0;0;2;20];
     %y0 = [27.8;0;1.22;2.62] %perfect KAM scenario
     y0 = [0.2,-0.2828,0,0];
     y0 = [pi;0;.5;.5]
     %y0 = [pi;pi;0.5;0]
     s=ode45(@doublependulum2,[0,T],y0,[],0,1,1,1,1);
     t = 0:0.01:T; x = deval(s,t)';
     %find energy for the inital conditions
     energy = hamiltonian(y0)
     %creates a plot with four figures of the variables plotted against each
     %other
     figure;
     subplot(222)
     plot(x(:,3),x(:,4));title('ang ular velocities')
     subplot(221)
     plot(x(:,1),x(:,2));title('ang les')
     subplot(223)
     plot(x(:,1),x(:,3));title('pos itions and momentum of 1')
     subplot(224)
     plot(x(:,2),x(:,4));title('pos itions and momentum of 1
    ')

    -----
    Dernière modification par yoda1234 ; 30/10/2011 à 11h56.

  2. #2
    invitea41a00b4

    Re : problème pour traduire un code Matlab en Mathematica

    aidez moi s'il vous plait c'est important, j'ai surtout du mal à traduire les boucles. Je vais apporter des précisions sur ce que je souhaite faire. Je travaille sur un double pendule. Les équations du mouvement sont données ici : http://gilbert.gastebois.pagesperso-...ule_double.htm. ( dans mon cas, m1=m2=l1=l2=1 ). J'ai déjà résolu les équations differentielles sur Mathematica et obtenu le graphe de la trajectoire et le portrait de phase, mais il me reste le plus important à programmer: la section de poincaré en 2D et 3D (1er code et 2ème code), et le graphe me donnant l'exposant de Lyapunov (3ème code). Merci d'avance

  3. #3
    lou_ibmix_xi

    Re : problème pour traduire un code Matlab en Mathematica

    Une reponse a cote de la plaque mais qui peut etre utile, tu peux utiliser "OCTAVE" qui est un clone gratuit de matlab, ainsi aucune traduction de scripts ne sera necessairre.

  4. #4
    invitea41a00b4

    Re : problème pour traduire un code Matlab en Mathematica

    Merci, je viens de le telecharger. Mais je ne comprends pas son fonctionnement :/ Comment fait-on ?

  5. A voir en vidéo sur Futura

Discussions similaires

  1. [matlab] programmation : peut-on générer un arbre des fonctions d'un code Matlab ?
    Par bratisla dans le forum Logiciel - Software - Open Source
    Réponses: 3
    Dernier message: 09/07/2015, 11h12
  2. [Matlab] Traduire la méthode de la Sécante
    Par invite78f958b1 dans le forum Programmation et langages, Algorithmique
    Réponses: 14
    Dernier message: 30/03/2014, 15h56
  3. conversion matlab vers mathematica
    Par invitea41a00b4 dans le forum Mathématiques du supérieur
    Réponses: 3
    Dernier message: 30/10/2011, 18h54
  4. problème avec ADPLL code matlab
    Par invite1d14401a dans le forum Logiciel - Software - Open Source
    Réponses: 0
    Dernier message: 02/11/2010, 10h44
  5. Mathematica CODE - Passer une fonction en paramètre
    Par Lévesque dans le forum Mathématiques du supérieur
    Réponses: 1
    Dernier message: 20/02/2006, 14h44