I have drawn the graph of nonlinear pendulum\'s equation. I need also in the sam

Question

I have drawn the graph of nonlinear pendulum's equation. I need also in the same plot the graph of the linear pendulum's equation to see the difference. Use matlab

the following code is for the nonlinear case

clear; clc;

L = 2; g=32.17;

F=@(t,Y) [ Y(2) ; ...

                   -(g/L)*sin(Y(1))];

t0=0; tf=2; n=1000; h=(tf-t0)/n;

alpha = [ pi/6 ; 0 ];

[t,w] = midpoint(t0,h,n,alpha,F);

[t',w']

plot(t,w(1,:));grid on

function [t,w] = midpoint(t0,h,n,alpha,F)

w(1:2,1)=alpha; t(1)=t0;

for i=1:n

    K1 = F(t(i),w(1:2, i));  

    K2 = F(t(i)+0.5*h,w(1:2, i)+0.5*h*K1);

    w(1:2, i+1)=w(1:2, i)+(K2)*h;

    t(i+1)=t(i)+h;

end

end

Explanation / Answer

function handles = plot_arrow(x1,y1,x2,y2,varargin)

%

% plot_arrow - plots an arrow to the current plot

% from: http://www.itl.nist.gov/div898/winds/wind_pressure/m-files/plot_arrow.m

%

% format: handles = plot_arrow( x1,y1,x2,y2 [,options...] )

%

% input: x1,y1 - starting point

% x2,y2 - end point

% options - come as pairs of "property","value" as defined for "line" and "patch"

% controls, see matlab help for listing of these properties.

% note that not all properties where added, one might add them at the end of this file.

%

% additional options are:

% 'headwidth': relative to complete arrow size, default value is 0.07

% 'headheight': relative to complete arrow size, default value is 0.15

% (encoded are maximal values if pixels, for the case that the arrow is very long)

%

% output: handles - handles of the graphical elements building the arrow

%

% Example: plot_arrow( -1,-1,15,12,'linewidth',2,'color',[0.5 0.5 0.5],'facecolor',[0.5 0.5 0.5] );

% plot_arrow( 0,0,5,4,'linewidth',2,'headwidth',0.25,'headheight',0.33 );

% plot_arrow; % will launch demo

% =============================================

% for debug - demo - can be erased

% =============================================

if (nargin==0)

figure;

axis;

set( gca,'nextplot','add' );

for x = 0:0.3:2*pi

color = [rand rand rand];

h = plot_arrow( 1,1,50*rand*cos(x),50*rand*sin(x),...

'color',color,'facecolor',color,'edgecolor',color );

set( h,'linewidth',2 );

end

hold off;

return

end

% =============================================

% end of for debug

% =============================================

% =============================================

% constants (can be edited)

% =============================================

alpha = 0.15; % head length

beta = 0.07; % head width

max_length = 22;

max_width = 10;

% =============================================

% check if head properties are given

% =============================================

% if ratio is always fixed, this section can be removed!

if ~isempty( varargin )

for c = 1:floor(length(varargin)/2)

try

switch lower(varargin{c*2-1})

% head properties - do nothing, since handled above already

case 'headheight',alpha = max( min( varargin{c*2},1 ),0.01 );

case 'headwidth', beta = max( min( varargin{c*2},1 ),0.01 );

end

catch

fprintf( 'unrecognized property or value for: %s ',varargin{c*2-1} );

end

end

end

% =============================================

% calculate the arrow head coordinates

% =============================================

den = x2 - x1 + eps; % make sure no devision by zero occurs

teta = atan( (y2-y1)/den ) + pi*(x2<x1) - pi/2; % angle of arrow

cs = cos(teta); % rotation matrix

ss = sin(teta);

R = [cs -ss;ss cs];

line_length = sqrt( (y2-y1)^2 + (x2-x1)^2 ); % sizes

head_length = min( line_length*alpha,max_length )*1.5;

head_width = min( line_length*beta,max_length )*1.5;

x0 = x2*cs + y2*ss; % build head coordinats

y0 = -x2*ss + y2*cs;

coords = R*[x0 x0+head_width/2 x0-head_width/2; y0 y0-head_length y0-head_length];

% =============================================

% plot arrow (= line + patch of a triangle)

% =============================================

h1 = plot( [x1,x2],[y1,y2],'k' );

h2 = patch( coords(1,:),coords(2,:),[0 0 0] );

% =============================================

% return handles

% =============================================

handles = [h1 h2];

% =============================================

% check if styling is required

% =============================================

% if no styling, this section can be removed!

if ~isempty( varargin )

for c = 1:floor(length(varargin)/2)

try

switch lower(varargin{c*2-1})

% only patch properties   

case 'edgecolor', set( h2,'EdgeColor',varargin{c*2} );

case 'facecolor', set( h2,'FaceColor',varargin{c*2} );

case 'facelighting',set( h2,'FaceLighting',varargin{c*2} );

case 'edgelighting',set( h2,'EdgeLighting',varargin{c*2} );

% only line properties   

case 'color' , set( h1,'Color',varargin{c*2} );

% shared properties   

case 'linestyle', set( handles,'LineStyle',varargin{c*2} );

case 'linewidth', set( handles,'LineWidth',varargin{c*2} );

case 'parent', set( handles,'parent',varargin{c*2} );

% head properties - do nothing, since handled above already

case 'headwidth',

case 'headheight',

end

catch

fprintf( 'unrecognized property or value for: %s ',varargin{c*2-1} );

end

end

end

Related Questions

Navigate

• Browse (All)
• Browse I
• Subjects

Integrity-first tutoring: explanations and feedback only — we do not complete graded work. Learn more.