Commit fdf7e180 authored by Jérémy Baudry's avatar Jérémy Baudry

nettoyage

parent 4e5b4414
clear all;
close all;
x=ncread('simulation1.nc','x_axis');
t=ncread('simulation1.nc','time');
om=ncread('simulation1.nc','omega');
spectre=ncread('simulation1.nc','Spectrum');
Dave=ncread('simulation1.nc','Dave');
Dmax=ncread('simulation1.nc','Dmax');
hice=ncread('simulation1.nc','Ice thickness');
c=ncread('simulation1.nc','Ice concentration');
FSD=ncread('simulation1.nc','Floe size distribution');
floesize=ncread('simulation1.nc','floe size');
% f=om/(2*pi);
% E=reshape(spectre(end,end,:),length(om),1);
% Ei=reshape(spectre(30,:,:),length(x),length(om));
% E1=reshape(spectre(1,1,:),length(om),1);
%
%
%
%
%
% figure(1)
% cmap=rand(length(t),3);
% w = waitforbuttonpress;
% for i=1:length(t)
%
% figure(1)
%
% sp=reshape(spectre(i,:,:),length(x),length(om));
% h=mesh(f,x,sp);
% axis([min(f) max(f) min(x) max(x) 0 max(E1)])
% xlabel('Frequency [s^{-1}]')
% zlabel('Energy')
%
% ylabel('x [km]')
% pause(0.1)
%
% end
%
%
%
%
% figure
%
% plot(x,hice,'color','r')
% hold on
% plot(x,c,'--b')
% grid on
%
% xlabel('x [km]')
%
% legend('ice thickness','ice concentration')
%
%
%
% figure
%
%
%
% plot(x,Dmax,'color','r')
% hold on
% plot(x,Dave,'--b')
% grid on
%
% xlabel('x [km]')
% ylabel('Floe size [m]')
%
% legend('Dmax','<D>')
\ No newline at end of file
This diff is collapsed.
subroutine break_horvat
use parameters
implicit none
integer::fl,iii,jjj,bin_new_floe
double precision :: rho=1025
double precision :: strain_crit=3e-5,new_floe
double precision, allocatable :: Amp(:),P_wa(:)
double precision, allocatable :: strain(:),P_f(:,:),coef(:)
allocate(Amp(nfreq))
allocate(P_wa(nfreq))
allocate(strain(nfreq))
allocate(P_f(nbcat,nfreq))
allocate(coef(nbcat))
P_f=0
!transform energy spectrum into amplitude spectrum
Amp(1:nfreq)=sqrt(2*E(n,i,1:nfreq)*domega/(rho*g))
!calcul the strain at each wavelength
strain(1:nfreq)=Amp(1:nfreq)*h*2*pi**2/wl(1:nfreq)**2
!find the probability of occurence of these amplitude using rayleigh
!distribution
P_wa(1:nfreq)=2*Amp(1:nfreq)*exp(-Amp(1:nfreq)**2/Hs**2)/Hs**2
!Compute breaking probability
fl=0
do iii=1,nedge-1
fl=fl+1
do jjj=1,nfreq
if (strain(jjj).gt.strain_crit .and. wl(jjj).lt.floe_cat(iii)) then
P_f(fl,jjj)=P_wa(jjj)
else
P_f(fl,jjj)=0
end if
end do
!coef(fl)=sum(P_f(fl,1:nfreq))
!P_f(fl,1:nfreq)=P_f(fl,1:nfreq)/coef(fl)
end do
! do the fragmentation
fl=0
do iii=1,nedge-1
fl=fl+1
do jjj=1,nfreq
if (strain(jjj).gt.strain_crit .and. wl(jjj).lt.floe_cat(iii)) then
new_floe=0.5*wl(jjj)
bin_new_floe=nedge-ceiling((new_floe-minfloe)/res)
FSD(n,i,bin_new_floe)=FSD(n,i,bin_new_floe)+P_f(fl,jjj)*FSD(n,i,fl)
FSD(n,i,fl)=FSD(n,i,fl)-P_f(fl,jjj)*FSD(n,i,fl)
end if
end do
end do
if(i.eq.2 .and. n.eq.2) then
do fl=1,nbcat
write(11,*)P_f(fl,1:nfreq)
end do
end if
end subroutine break_horvat
subroutine graph_dislin
use parameters
use dislin
IMPLICIT NONE
double precision, allocatable :: axe(:),Ep(:,:)
real ::par1,par2,par3,par4,par5,par6,par7,par8,par9,par10,par11,par12
allocate(axe(nbin))
allocate(Ep(nbin,nfreq))
axe(1)=0
do i=1,nbin-1
axe(i+1)=i*dx/1000
end do
!_____________________________________maximum floe size________________
!CALL SETPAG('DA4P')
!CALL METAFL('CONS')
!CALL SCRMOD('REVERS')
!CALL DISINI()
!CALL PAGERA()
!CALL SETRGB (0.5,0.2,0.4)
!CALL GRAF (real(axe(1),4),real(axe(nbin),4),real(0,4),real(10*dx/1000,4),real(0,4),real(maxval(Dmax),4),real(0,4),real(50,4))
!call NAME('x axis [km]','X')
!call NAME('Dmax [m]','Y')
!CALL TITLIN ('Maximum Floe size in the MIZ',1)
!call TITLE
! call COLOR('GRAY')
! call GRID(10,5)
! call COLOR('RED')
! CALL CURVE (real(axe,4),real(Dmax,4),nbin)
! call ENDGRF
!call disfin
!________________________SPECTRUM_______________________________________
do i=1,nbin
do ii=1,nfreq
Ep(i,ii)=real(E(15,i,ii),4)
end do
end do
par5=real(omega(1),4)
par6=real(omega(nfreq),4)
par7=real(omega(1),4)
par8=real(1,4)
par1=real(axe(1),4)
par2=real(axe(nbin),4)
par3=real(0,4)
par4=real(10*dx/1000,4)
par9=real(0.,4)
par10=real(maxval(Ei),4)
par11=real(0.,4)
par12=real(maxval(Ei)/10,4)
!write(*,*)par1,par2,par3,par4,par5,par6,par7,par8,par9,par10,par11,par12
CALL SETPAG('DA4P')
CALL METAFL('CONS')
!CALL SCRMOD('REVERS')
CALL DISINI()
CALL PAGERA()
CALL VIEW3D(-5.,-5.,4.,'ABS')
do n=1,nsteps
do i=1,nbin
do ii=1,nfreq
Ep(i,ii)=real(E(n,i,ii),4)
end do
end do
call ERASE
call COLOR('WHITE')
CALL GRAF3D (par1,par2,par3,par4,par5,par6,par7,par8,par9,par10,par11,par12)
CALL NAME('X-axis [km]','X')
CALL NAME('Frequency','Y')
CALL NAME('Energy','Z')
call COLOR('BLUE')
CALL SURMAT (real(Ep,4),nbin,nfreq,1,1)
call ENDGRF
call system('sleep 0.1')
end do
call disfin
end subroutine graph_dislin
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment