subroutine attenuation use parameters implicit none double precision :: q11,q12,q13,q14,q15,q21,q22,q23,q24,q25 double precision :: q31,q32,q33,q34,q35 double precision, allocatable ::p1(:),p2(:),p3(:),att(:) allocate(p1(nfreq)) allocate(p2(nfreq)) allocate(p3(nfreq)) allocate(att(nfreq)) q11 = -0.00000777 q12 = 0.00032080 q13 = -0.00437542 q14 = 0.02047559 q15 = -0.01356537 q21 = 0.00003635 q22 = -0.00153484 q23 = 0.02121709 q24 = -0.09289399 q25 = -0.03693082 q31 = -0.00004509 q32 = 0.00214484 q33 = -0.03663425 q34 = 0.26065369 q35 = -0.62474085 p1 = q11*T**4 + q12*T**3 + q13*T**2 + q14*T + q15 p2 = q21*T**4 + q22*T**3 + q23*T**2 + q24*T + q25 p3 = q31*T**4 + q32*T**3 + q33*T**2 + q34*T + q35 alpha(i,1:nfreq)=-1*(p1*h(i)**2 + p2*h(i) + p3) do j=1,nfreq if (alpha(i,j).lt.0)then alpha(i,j)=0d0 end if end do att=C_ice(i)*alpha(i,1:nfreq)/(Dave(i)+3e-14) S_ice(i,1:nfreq)=-att*E(n,i,1:nfreq) E(n,i,1:nfreq)=E(n,i,1:nfreq)*exp(-att*Cg*dt) end subroutine attenuation