Skip to main content

Posts

Second Harmonic Generation N=1:21

gnuplot> set xrange [-180:180] gnuplot> set yrange [-180:180] splot sin(cos(x*pi/180))*sin(cos(x* pi/180))/(cos(x*pi/180)*cos(x* pi/180))*sin(cos(y*pi/180))* sin(cos(y*pi/180))/(cos(y*pi/ 180)*cos(y*pi/180))+sin(cos(( x+1)*pi/180))*sin(cos((x+1)* pi/180))/(cos((x+1)*pi/180)* cos((x+1)*pi/180))*sin(cos((y+ 1)*pi/180))*sin(cos((y+1)*pi/ 180))/(cos((y+1)*pi/180)*cos(( y+1)*pi/180))+sin(cos((x+2)* pi/180))*sin(cos((x+2)*pi/180) )/(cos((x+2)*pi/180)*cos((x+2) *pi/180))*sin(cos((y+2)*pi/ 180))*sin(cos((y+2)*pi/180))/( cos((y+2)*pi/180)*cos((y+2)* pi/180))+sin(cos((x+3)*pi/180) )*sin(cos((x+3)*pi/180))/(cos( (x+3)*pi/180)*cos((x+3)*pi/ 180))*sin(cos((y+3)*pi/180))* sin(cos((y+3)*pi/180))/(cos(( y+3)*pi/180)*cos((y+3)*pi/180) )+sin(cos((x+4)*pi/180))*sin( cos((x+4)*pi/180))/(cos((x+4)* pi/180)*cos((x+4)*pi/180))* sin(cos((y+4)*pi/180))*sin( cos((y+4)*pi/180))/(cos((y+4)* pi/180)*cos((y+4)*pi/180)) +sin(cos((x+5)*pi/180))*sin( cos((x+5)*pi/180))/(cos((x+5)* pi/180)*cos((
Recent posts

Second Harmonic Generation

gnuplot> set xrange [-180:180] gnuplot> set yrange [-180:180] gnuplot> set pm3d gnuplot> set hidden3d  gnuplot> set title 'SHG' gnuplot> splot sin(cos(x*pi/180))*sin(cos(x* pi/180))/(cos(x*pi/180)*cos(x* pi/180))*sin(cos(y*pi/180))* sin(cos(y*pi/180))/(cos(y*pi/ 180)*cos(y*pi/180)) title 'N=1 '

FREE FALL

## A function, free fall, that takes in h (in metres), and ## returns the final velocity of the ball at the ## time step just before it touches the ground and makes ## a plot of velocity versus time. function freefall(h) ## constants and initializations v=0; ## initial velocity(at rest) [m/sec] y=h; ## initial altitude of the ball from the ground [m] t=0; ## initial time[sec ] B1=0.05; ## Coeff of the term prop to v [kg/sec] B2=6E-4; ## Coeff of the term prop to v^2 [kg/m] m=0.25; ## Mass of the ball [kg] g=9.8; ## Gravitational acceleration [m/sec^2] dt=0.1; ## Time increment [sec] n=1; ## Initialize the loop index[dimensionless] ## Run the loop or iterate until the vertical component is ## smaller than zero. while (y(n)>0); ## decrease the altitude in each step y=[y;y(n)+v(n)*dt]; ## and accumlate the results in the array ## of the vertical displacement. ## increase the time in each step t=[t;t(n)+dt]; ## and accumlate the results in the array ## of the tim

INTEGRAL (TRAPEZOIDAL RULE)

## This m-file computes the integral of a function f over the ## equally spaced points x using the trapezoidal rule. function I = trapezoid(x,f) h = x(2)-x(1); N = length(f); I = 0; for n=2:N-1, I = I + h*f(n); endfor I = I + h/2*(f(1) + f(N)); endfunction

TOTAL ISING ENERGY

function E=Ising_Energy(lattice,i,j,rneigh) H=2; ## applied magnetic field times permeability[joule] E=0; ## initial configuration energy[joule] J0=1 ## coupling constant of function J, normaly in joule. alpha=2 ## exponential constant of function J. L=length(lattice); ## number of particles or spins in lattice for r=1:rneigh NN1=mod(j+r,L); NN1+=(NN1==0)*L; NN2=mod(j-r,L); NN2+=(NN2==0)*L; NN3=mod(i+r,L); NN3+=(NN3==0)*L; NN4=mod(i-r,L); NN4+=(NN4==0)*L; J=coupling(rneigh,J0,alpha); dE=J*(lattice(i,NN1)+lattice(i,NN2)+ lattice(NN3,j)+lattice(NN4,j))-H*(sum((rand(L)>0.5)*2-1)); E+=-1*dE endfor endfunction

NONUNFORM SOUND WAVE GENERATOR

## a script for a sound wave source made up of two different masses ##seperated in the middle. Let the rigth part is lighter than the ##left one that ensures that the ligther is faster, and so has ##greater r since r=cdt/dx which is dimensionless(ratio). dx=1e-2; ## Increment in the horizontal distance (m) c1=300; ## speed of the left one (m/s) c2=500; ## Speed of the rigth one (m/s) dt=dx/max(c1,c2); ## Increment in time (s). The max c, the min dt. r1=c1*dt/dx; ## Ratio r for the left one r2=c2*dt/dx; ## Ration r for the right one x=-1:dx:1; ## same as the string_fixed l=length(x); x0=0.5; k=1e2; ## Inital profile(return to guassuian distibution) y=initial_profile(x,x0,k); axis([-1.05,1.05,-1.1,1.1]) plot(x,y,'r;;',[0 0],[-1.1 1.1],'r-;;') pause ## Boundary conditions while t=dtn=0 ynow=y; yprev=y; N=100; for n=1:N ynext=propagate_two_parts(ynow,yprev,r1,r2); ## calling the prepared ## function plot(x,ynext,';;',[0 0],[-1.1 1.1],'r