upm.m

是有分布傅立叶算法解薛定谔方程的一个程序（下面我就不乱说说了请认真阅读您的文件包然后写出其具体功能(至少要20个字)。尽量不要让站长把时间都花费在为您修正说明上。压缩包解压时不能有密码。系统会自动删除

function u1 = UPM(u0,dt,dz,nz,alpha,betap,gamma,tol);


% This function solves the nonlinear Schrodinger equation for
% pulse propagation in an optical fiber using the split-step
% Fourier method described in: 
% 
%  A. A. Rieznik,  T. Tolisano,  F. A. Callegari,  D. F. Grosz, and H. L. Fragnito, 
% "Uncertainty relation for the optimization of optical-fiber transmission
% systems simulations," Opt. Express 13, 3822-3834 (2005).
% http://www.opticsexpress.org/abstract.cfm?URI=OPEX-13-10-3822
%
% USAGE
%
% u1 = ssprop(u0,dt,dz,nz,alpha,betap,gamma);
% u1 = ssprop(u0,dt,dz,nz,alpha,betap,gamma,tol);
%
% INPUT
%
% u0 - starting field amplitude (vector)
% dt - time step
% dz - propagation stepsize
% nz - number of steps to take, ie, ztotal = dz*nz
% alpha - power loss coefficient, ie, P=P0*exp(-alpha*z)
% betap - dispersion polynomial coefs, [beta_0 ... beta_m]
% gamma - nonlinearity coefficient
% tol - convergence tolerance (default = 1e-5)
%
% OUTPUT
%
% u1 - field at the output
% 
% NOTES  The dimensions of the input and output quantities can
% be anything, as long as they are self consistent.  E.g., if
% |u|^2 has dimensions of Watts and dz has dimensions of
% meters, then gamma should be specified in W^-1*m^-1.
% Similarly, if dt is given in picoseconds, and dz is given in
% meters, then beta(n) should have dimensions of ps^(n-1)/m.

if (nargin<8)
  tol = 1e-3;
end

nt = length(u0);
w = 2*pi*[(0:nt/2-1),(-nt/2:-1)]'/(dt*nt);

linearoperator = -alpha/2;
for ii = 0:length(betap)-1;
  linearoperator = linearoperator - j*betap(ii+1)*(w).^ii/factorial(ii);
end

lineraoperatorwithoutloss = linearoperator + alpha/2;

u1 = u0;
ufft = fft(u0);

fiberlength = nz*dz;
propagedlength =0;

fprintf(1, '\nSimulation running...      ');


while propagedlength < fiberlength,
    
  Et = sum(abs(u1).^2);
  meanN = sum(gamma*(abs(u1).^4)) / Et;
  aux = (gamma*(abs(u1).^2)-meanN).^2;
  deltaN = sqrt(sum( aux.*abs(u1).^2 ) / Et);
  
  Ef = sum(abs(ufft).^2);
  meanD = sum( j*lineraoperatorwithoutloss.*(abs(ufft).^2) ) / Ef;
  aux = ( j*lineraoperatorwithoutloss - meanD ).^2;
  deltaD = sqrt(sum( aux.*abs(ufft).^2 ) / Ef);
  
   
  dz = (tol^(1/2))*sqrt(1 / (deltaD*deltaN));
  
  if (dz + propagedlength) > fiberlength,
      dz = fiberlength - propagedlength;
  end
  
  halfstep = exp(linearoperator*dz/2);
  uhalf = ifft(halfstep.*ufft);
  
  u1 = uhalf .* exp(-j*gamma*(abs(uhalf).^2 )*dz);
  ufft = halfstep.*fft(u1);

  propagedlength = propagedlength + dz;
  fprintf(1, '\b\b\b\b\b\b%5.1f%%', propagedlength * 100.0 /fiberlength );
end

ux = ifft(ufft);
u1=ux;