瑞雷面波反演程序

瑞雷面波反演，很好的程序，物探人员可以学习一下，使用的是基因算法。反演结果准确，且程序比较清楚。

function [zrdr] = disp_stress（freqvrthkdnscvscvpFz）

% This function calulates the displacement-stress vectors （i.e. the eigenfunctions）
% corresponding to the phase velocities （i.e. wavenumbers） contained in the vr
% matrix

% Copyright 1999 by Glenn J. Rix and Carlo G. Lai

% This program is free software; you can redistribute it and/or
% modify it under the terms of the GNU General Public License
% as published by the Free Software Foundation; either version 2
% of the License or （at your option） any later version.
%
% This program is distributed in the hope that it will be useful
% but WITHOUT ANY WARRANTY; without even the implied warranty of
% MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
% GNU General Public License for more details.
%
% You should have received a copy of the GNU General Public License
% along with this program; if not write to the Free Software
% Foundation Inc. 59 Temple Place - Suite 330 Boston MA  02111-1307 USA.

% Establish global parameters
global NUMPOINTS LAMBDA MAXROOT

% Calculate the vector of circular frequencies
om = 2*pi*freq;

% Determine the number of layers not including the half space
N = length（thk）;

% Calulate the maximum depth for determining the displacement-stress vectors
lambda_max = LAMBDA*（real（vr（:1））.^2 + imag（vr（:1））.^2）./real（vr（:1））./freq;

% Initiate the depth and displacement-stress vectors and their numerical derivatives
z = zeros（NUMPOINTSlength（freq））;
r = zeros（length（freq）MAXROOTNUMPOINTS4）;
dr = zeros（length（freq）MAXROOTNUMPOINTS2）;

% Loop through the frequencies
for j = 1:length（freq）
   
   % Create a vector of depths
   z（:j） = linspace（0lambda_max（j）NUMPOINTS）‘;
   
   % Loop through the modes at each frequency
   index1 = find（vr（j:））;
   for m = 1:length（index1）
      
      % Calculate the wavenumber and load vector
      k = om（j）/vr（jindex1（m））;
      delqz = [0 ; k*Fz/（2*pi）];
      
      % Check to see if the phase velocity is equal to the shear wave velocity
		% or compression wave velocity of one of the layers
      epsilon = 0.0001;
      while any（abs（om（j）/k-cvs）   		k = k * （1+epsilon）;
		end   
      
      % Calculate the PSV element matrices for each layer and generalized R/T matrices
		[e11e12e21e22dumunusnup] = psv（thkdnscvpcvsom（j）k）;
		[tdturdru] = modrt（e11e12e21e22du）;
		[TdRd] = genrt（tdturdru）;
      
      % Initialize the Cd and Cu matrices
      cd = zeros（21N+1）;
      cu = zeros（21N+1）;
      
      % Calculate Cd for the first layer
      [lamdlamu] = updown（thkcvpcvsom（j）k01）;
      cd（::1） = （e21（::1） + e22（::1）*lamu*Rd（::1））\delqz;
      
      % Calculate Cd and Cu for the remaining layers
      for n = 1:N
         cu（::n） = Rd（::n）*cd（::n）;
         cd（::n+1） = Td（::n）*cd（::n）;
      end
      
      % Loop through the ve

 属性            大小     日期    时间   名称
----------- ---------  ---------- -----  ----

     文件       3217  2005-08-11 16:26  GA-based_Inversion_Example1-1\BS2RV.M

     文件       1835  2005-08-11 16:26  GA-based_Inversion_Example1-1\CONTENTS.M

     文件       1168  2005-08-11 16:26  GA-based_Inversion_Example1-1\CRTbase.M

     文件       2273  2006-12-25 13:54  GA-based_Inversion_Example1-1\CRTBP.M

     文件       2091  2005-08-11 16:26  GA-based_Inversion_Example1-1\CRTRP.M

     文件       3884  1999-12-23 16:20  GA-based_Inversion_Example1-1\disp_stress.m

     文件       4193  2000-01-03 08:56  GA-based_Inversion_Example1-1\effective.m

     文件       9684  2008-07-10 12:40  GA-based_Inversion_Example1-1\Example_1.asv

     文件       9684  2008-07-03 15:56  GA-based_Inversion_Example1-1\Example_1.m

     文件       7448  2008-06-22 19:44  GA-based_Inversion_Example1-1\Figure No.1.fig

     文件      21992  2008-06-22 19:44  GA-based_Inversion_Example1-1\Figure No.2.fig

     文件      20976  2008-06-22 19:44  GA-based_Inversion_Example1-1\Figure No.3.fig

     文件       1270  1999-11-11 08:16  GA-based_Inversion_Example1-1\genrt.m

     文件       2429  1999-12-23 13:40  GA-based_Inversion_Example1-1\green.m

     文件       1475  1999-11-11 08:17  GA-based_Inversion_Example1-1\homogeneous.m

     文件       7205  2005-08-11 16:26  GA-based_Inversion_Example1-1\MIGRATE.M

     文件       1622  2006-12-28 13:26  GA-based_Inversion_Example1-1\Mlayer_NumEx1.m

     文件       2698  1999-11-11 08:17  GA-based_Inversion_Example1-1\modal.m

     文件       1548  1999-11-11 08:17  GA-based_Inversion_Example1-1\modrt.m

     文件       4019  2005-08-11 16:26  GA-based_Inversion_Example1-1\MPGA.M

     文件       1609  2005-08-11 16:26  GA-based_Inversion_Example1-1\MUT.M

     文件       3437  2005-08-11 16:26  GA-based_Inversion_Example1-1\MUTATE.M

     文件       4943  2005-08-11 16:26  GA-based_Inversion_Example1-1\MUTBGA.M

     文件        393  2005-10-07 13:16  GA-based_Inversion_Example1-1\NumExample1.txt

     文件       2629  2005-08-11 16:26  GA-based_Inversion_Example1-1\OBJFUN1.M

     文件       2559  2005-08-11 16:26  GA-based_Inversion_Example1-1\OBJHARV.M

     文件       4187  1999-12-24 16:35  GA-based_Inversion_Example1-1\partial.m

     文件       2896  1999-11-11 08:18  GA-based_Inversion_Example1-1\psv.m

     文件       4709  2005-08-11 16:26  GA-based_Inversion_Example1-1\RANKING.M

     文件       7593  2006-12-25 13:50  GA-based_Inversion_Example1-1\rayleigh.m

............此处省略26个文件信息