Skip to main content

REALISTIC AND NON-REALISTIC STRINGS


## A 'realistic', 'non-elastic' string, which responses to any
## bending and has stifness. This script takes in the previous
## and the present profiles and iterates to find
## the profile in the next time step. The ratio 'r' is not 1
## like in the 'non-realistic' string since the speed of the wave
## always less than the speed of the string it should be less than 1
## for best and most stable solution
##constants
dx=1e-2 ## Spatial increment (m)
L=2 ## Length of the string (m)
M=L/dx ## Dimensionless partition
E=1e-4 ## Dimensionless stiffnes
function ynext=propagate_stiff(ynow,yprev,r)
## Quick and dirty way to fix boundary conditions -- for each step
## they are the same as the previous step.
ynext=ynow;
ynow(1)=ynow(2)=0;
##Entering the loop
for i=3:length(ynow)-1
## boundary condition
ynow(length(ynow)-1)=ynow(length(ynow)-2)=0;
## Divide the ynext with many terms into three parts for easiness
ynext(i)=(2−(2*r^2)−(6*E*(r^2)*(M^2)))*ynow(i)−yprev(i);
ynext(i)=ynext(i)+(r^2)*(1+4*E*(M^2))*(ynow(i+1)+ynow(i−1));
ynext(i)=ynext(i)-E*(r^2)*(M^2)*(ynow(i+2)+ynow(i−2));
endfor
endfunction

Comments

Popular posts from this blog

NEWTON-RAPSON METHOD FOR HEAT FLOW

##Constants and initializations
a=5.67E-8; ## Stefan-Boltzman constant[Watt/meter^2Kelvin^4]
e=0.8; ## Rod surface emissivity [Dimensionless]
h=20; ## Heat transfer coefficient of air flow [W/m^2-K]
Tinf=Ts=25; ## Temperature of air and the walls of the close[Celcius]
D=0.1; ## Diameter of the rod[meter]
I2R=100; ## Electric power dissipated in rod (Ohmic Heat)[W]
T=[]; ## Temperature of the rod[*C]
T(1)=25; ## Initial guess of the temperature of the rod[*C]
Q=[]; ## Heat function [W]
Qp=[]; ## First derivative of Q wrt T [W/C*].
for i=1:100
Q(i)=pi*D*(h*(T(i)-Tinf)+e*a*(T(i)^4-Ts^4))-I2R;
Qp(i)=pi*D*(h+4*e*a*T(i)^3);
T(i+1)=T(i)-Q(i)/Qp(i); ## Newton-Rapson Method
endfor
printf('The steady state temperature is %f\n',T(i+1))
save -text HeatFlowTemp.dat
## The plot
t=1:100; ##temperature
for n=1:100
H(n)=pi*D*(h*(t(n)-Tinf)+e*a*(t(n)^4-Ts^4))-I2R;
endfor
plot(t,H)
xlabel('T(Celcius)');
ylabel('Q(Watt)');
legend('Q(T)');
title('Heat flow vs Temperatu…