Newsgroups: bit.listserv.sas-l
Date:         Wed, 23 Dec 1992 08:30:03 CET
Reply-To:     Arnold Schick <SCHICK@DMRHRZ11.BITNET>
From:         Arnold Schick <SCHICK@DMRHRZ11.BITNET>
Subject:      differentation in case of 3D

Hallo SAS-Lers,

since I had sent a %macro to differentiate a function y=f(x), now I would
like to present a %macro to solve partial differentation of functions with
two variables as z=f(x,y). It's also a present to _ALL_ for christmas and
a big Thanks for your helps in the past.

The theory of the differentation is like this:

         Y                              parameter-formula:
          ^           / f(x)
   y(i+2) |          +                  y'   = (y    - y )/(x   - x ) + ERR
          |        /                     i+1     i+1    i    i+1   i
   y(i+1) |      +
          |     /
     y(i) |    +                        with y'    = 2*y' - y'
          |                                    i=1      3    2
          +------------------------>
              x(i)   x(i+2)        X    from i=1 to n points.
                 x(i+1)                 With ERR function equal zero by
                                        small delta X = x(i+1) - x(i).

 And in case of z=f(x,y) is to differentiate with this chema partial in steps,
 with varying of the dependent variable. The total derivation is the sum of
 differential to X and to Y.

The %macro DIFFER3D handles a differentation for z=f(x,z) with 3 SAS
data set variables. The %macro call is:

       %differ3d (data_in, data_out, x, y, z, m, nosort);

 x,y,z names the SAS data set variables within DATA_IN, which contains
 the parameter of function z=f(x,y). If data set DATA_IN is not present
 or missing, data set _LAST_ is inuse.

 the data set DATA_OUT includes the orginal x,y,z values and the z_diff
 variable with them values. If data set DATA_OUT is not present, halts
 the %macro.

 the factor m is a level of the source function: if m=0 the source function
 is not modified, if m=negative the source function becomes a fall
 (inclination) to the positive dependent variable (i.e. X) direction, and
 if m=positive rise it in this direction. Factor m is to use for a correction
 of the function z=f(x,y) like z=m*dependent_variable by differentation.
 If m is not present or missing, m=0 is inuse.

 the SORT option tells the %macro, if the input data set is not sorted
 by  dependent variables, and it is to sort through the %macro, NOSORT
 depresses it. If this option is not present, NOSORT is inuse.

In cases of needed inter-values between two points i and i+1 is to fit
with the G3GRID procedure before the %macro was invoked.


So long,

merry Christmas and a Happy N.Y. 1993

Arnold Schick  University of Marburg/Germany
----------------------------------------------------------------------
Here is this %MACRO for your use with 3 SAS data set examples (3rd example
is inuse and the others are enclosed in a comment):

%macro differ3d (data, out, x, y, z, m, sorts); *developed by Arnold Schick;
options nosource nostimer nonotes nosymbolgen;  *ACC, University of Marburg;

  %if &data  =  or &data  = . %then %let data = _LAST_ ;
  %if &out   =  or &out   = . %then %goto quit;
  data _null_;
    b = symget ('m') / 1;
    if b = . then b = 0;
    call symput('m',b);
  run;

  %if &sorts = sort or &sorts = . or &sorts =  %then %do;
      proc sort data = &data  out = &out;
        by &x &y;
      run;
    %end;
      %else %do;
         data &out;
           set &data;
         run;
      %end;

  data _null_ ;
    set &out;
    if _N_ = 1 then call symput('y_1', &y);
    if _N_ = 2 then call symput('y_2', &y);
    if _N_ = 3 then do; call symput('y_3', &y); stop; end;
  run;

  data anfang(keep= &x dz_start);
     set &out ;
     x1 = lag( &x );
     y1 = lag( &y );
     z1 = lag( &z ) + &m * &y;
     if y1 = &y then delete;
     if &y_2 = &y or &y_3 = &y then do;
         dz   = ( &z + &m * &y - z1 )/( &y - y1 );
         dz_b = lag(dz);
         if &y_3 = &y then do;
            dz_start = 2*dz_b - dz;
            output anfang;
         end;
     end;
  run;

  data &out;
    set &out;
    merge &out anfang;
       by &x;
    length z_diff 8 ;
    y_previo = lag( &y );
    z_previo = lag( &z ) + &m * &y;
    if &y = y_previo then delete;
    if &y_1 = &y then z_diff = dz_start;
                 else if &m ^= 0 then z_diff = -(z_previo - &z + &y * &m)
                                               /( &y - y_previo);
                                 else z_diff = -(z_previo - &z)
                                               /( &y - y_previo);
    drop y_previo z_previo dz_start;
  run;

  %goto final ;
  %quit : %put HALT: please define OUT data set (2nd parameter) ;
  %final : ;

options source stimer notes;
%mend differ3d;

/*
* SAS data set example one: surface z=exp(-x*y)/(2-cos(x*y))*sin(y-x);
data eins;
  do x=0 to 1 by 0.05;
     do y=0 to 2 by 0.1;
        z=exp(-x*y)/(2-cos(x*y))*sin(y-x);
        output;
     end;
  end;
run;

* SAS data set example two: Cowboy Hat  z=sin(sqrt(x**2 + y**2));
data eins;
  do x=-5 to 5 by 0.25;
     do y=-5 to 5 by 0.25;
        z = sin(sqrt(x*x + y*y));
        output;
     end;
  end;
run;
*/
* SAS data set example three: reconstruction of the Cowboy Hat;
data eins;
 do i=1 to 30;
  x=10*ranuni(33)-5;
  y=10*ranuni(35)-5;
  z=sin(sqrt(x*x + y*y));
  output;
 end;
run;
proc g3grid data=eins out=eins;
  grid x*y=z / spline
               axis1=-5 to 5 by 0.25
               axis2=-5 to 5 by 0.25;
run;
proc sort data=eins out=eins;
  by x y;
run;

%differ3d (eins, zwei,  x, y, z     ,-0.0025,nosort);
%differ3d (eins, drei,  y, x, z     ,.,sort);
%differ3d (drei, sechs, x, y, z_diff,0,sort);
%differ3d (zwei, fuenf, y, x, z_diff,,sort);

* the following part shows the results;
title 'orginal function';
proc g3d data=eins;
  plot y*x=z;
run;
title 'partial derivative to X';
proc g3d data=zwei;
  plot y*x=z_diff ;
run;
title 'partial derivative to Y';
proc g3d data=drei;
  plot y*x=z_diff;
run;
title 'partial derivative to XY';
proc g3d data=fuenf;
  plot y*x=z_diff;
run;
title 'partial derivative to YX';
proc g3d data=sechs;
  plot y*x=z_diff;
run;

*---------this part counts the total derivative and shows the result;
proc sort data=drei out=drei; by x y; run;
data vier;
  set drei (rename=(z_diff=zz_diff));
  merge drei zwei;
    by x y;
  z_diff = zz_diff+z_diff;
run;
title 'total derivative to X and Y';
proc g3d data=vier;
  plot y*x=z_diff;
run;
Comments: Gated by NETNEWS@AUVM.AMERICAN.EDU