Octave Encyclopedia

newmark_int

Statistics + Control & Dynamics 중심의 Octave 학습 예제

ex-04/newmark_int.m

코드 인덱스로 돌아가기

카테고리

Exercise 04

학습 소스 그룹

코드 길이

80

lines

작성자

-

날짜 정보 없음

패키지

none

pkg load 기준

함수 시그니처

function [u, udot, u2dot] = newmark_int(t, p, u0, udot0, m, k, xi, varargin)

전체 코드

전체 코드를 복사해서 Octave에서 바로 실행할 수 있습니다.

function [u, udot, u2dot] = newmark_int(t,p,u0,udot0,m,k,xi,varargin)
%Newmark's Direct Integration Method
%--------------------------------------------------------------------------
% Integrates a 1-DOF system with mass "m", spring stiffness "k" and damping
% coeffiecient "xi", when subjected to an external load P(t).
% Returns the displacement, velocity and acceleration of the system with
% respect to an inertial frame of reference.
%
% SYNTAX
%       [u, udot, u2dot] = newmark_int(t,p,u0,udot0,m,k,xi,varargin)
%
% INPUT
%       [t] :       Time Vector             [n,1]
%       [p] :       Externally Applied Load [n,1]
%       [u0]:       Initial Position        [1,1]
%       [udot0]:    Initial Velocity        [1,1]
%       [m]:        System Mass             [1,1]
%       [k]:        System Stiffness        [1,1]
%       [xi]:       System Damping          [1,1]
%       [varargin]: Options
%
% OUTPUT
%       [u]:        Displacemente Response   [n,1]
%       [udot]:     Velocity                 [n,1]
%       [u2dot]:    Acceleration             [n,1]
%
%  N = number of time steps
%
% The options include changing the value of the "gamma" and "beta"
% coefficient which appear in the formulation of the method. By default
% these values are set to gamma = 1/2 and alpha = 1/4.
%
% EXAMPLE
% To change nemark's coefficients, say to gamma = 1/3 and alpha = 1/5, 
% the syntax is:
%       [u, udot, u2dot] = newmark_int(t,p,u0,udot0,m,k,xi, 1/3, 1/5)  
%
%==========================================================================
%                     2004 By: Jose Antonio Abell Mena  (ja.abell@gmail.com)
if nargin == 7
%     disp('Using default values:');
%     disp('    gamma = 1/2');
%     disp('    beta  = 1/4');
    gam = 1/2;
    beta = 1/4;
else
    if nargin == 9
        beta = varargin{1};
        gam = varargin{2};
        
    else
        error('Incorrect number of imput arguments');
    end
end
wn = sqrt(k/m);
% wd = 2*xi*wn;

dt = t(2) - t(1);
c = 2*xi*wn*m;

kgor = k + gam/(beta*dt)*c + m/(beta*dt^2);
a = m/(beta*dt) + gam*c/beta;
b = 0.5*m/beta + dt*(0.5*gam/beta - 1)*c;

dp = diff(p);
u = zeros(length(t),1); udot = u; u2dot = u;
u(1) = u0;
udot(1) = udot0;
u2dot(1) = 1/m*(p(1)-k*u0-c*udot0);

for i = 1:(length(t)-1)
    deltaP = dp(i) + a*udot(i) + b*u2dot(i);
    du_i = deltaP/kgor;
    dudot_i = gam/(beta*dt)*du_i - gam/beta*udot(i) + dt*(1-0.5*gam/beta)*u2dot(i);
    du2dot_i = 1/(beta*dt^2)*du_i - 1/(beta*dt)*udot(i) - 0.5/beta*u2dot(i);
    u(i+1) = du_i + u(i);
    udot(i+1) = dudot_i + udot(i);
    u2dot(i+1) = du2dot_i + u2dot(i);
end

코드 해설

목적

Statistics + Control & Dynamics 중심의 Octave 학습 예제

입력

  • 파라미터: t
  • 파라미터: p
  • 파라미터: u0
  • 파라미터: udot0
  • 파라미터: m
  • 파라미터: k
  • 파라미터: xi
  • 파라미터: varargin

출력

  • 반환값: u
  • 반환값: udot
  • 반환값: u2dot
  • 콘솔 텍스트 출력

실행 흐름

  1. Integrates a 1-DOF system with mass "m", spring stiffness "k" and damping

핵심 함수/주제

u2dotudotdisplengthdiffdperrorMena

실습 과제

  • 같은 연산을 내장 함수와 사용자 함수 두 방식으로 계산해 오차를 비교해보세요.
  • 질량/감쇠/강성 또는 전달함수 계수를 바꿔 응답 변화를 확인해보세요.
  • 핵심 함수 u2dot의 인자를 한 가지 바꿔 결과 변화를 기록해보세요.

학습 팁

  • 입력 파일 경로가 현재 작업 디렉터리 기준인지 먼저 확인하세요.

같은 카테고리의 다른 코드