function [ suborder_xyzz, suborder_w ] = keast_subrule ( rule, suborder_num )

%*****************************************************************************80
%
%% KEAST_SUBRULE returns a compressed Keast rule.
%
%  Licensing:
%
%    This code is distributed under the GNU LGPL license.
%
%  Modified:
%
%    25 June 2007
%
%  Author:
%
%    John Burkardt
%
%  Reference:
%
%    Patrick Keast,
%    Moderate Degree Tetrahedral Quadrature Formulas,
%    Computer Methods in Applied Mechanics and Engineering,
%    Volume 55, Number 3, May 1986, pages 339-348.
%
%  Parameters:
%
%    Input, integer RULE, the index of the rule.
%
%    Input, integer SUBORDER_NUM, the number of suborders of the rule.
%
%    Output, real SUBORDER_XYZZ(4,SUBORDER_NUM),
%    the barycentric coordinates of the abscissas.
%
%    Output, real SUBORDER_W(SUBORDER_NUM), the
%    suborder weights.
%
  if ( rule == 1 )

    suborder_xyzz(1:4,1:suborder_num) = [ ...
      0.250000000000000000,  0.250000000000000000  ...
      0.250000000000000000,  0.250000000000000000  ...
    ]';

    suborder_w(1:suborder_num) = [ ...
      0.166666666666666667 ...
    ];

  elseif ( rule == 2 )

    suborder_xyzz(1:4,1:suborder_num) = [ ...
      0.585410196624968500,  0.138196601125010500  ...
      0.138196601125010500,  0.138196601125010500  ...
    ]';

    suborder_w(1:suborder_num) = [ ...
      0.0416666666666666667 ...
    ];

  elseif ( rule == 3 )

    suborder_xyzz(1:4,1:suborder_num) = [ ...
      0.250000000000000000,  0.250000000000000000  ...
      0.250000000000000000,  0.250000000000000000; ...
      0.500000000000000000,  0.166666666666666667  ...
      0.166666666666666667,  0.166666666666666667; ...
    ]';

    suborder_w(1:suborder_num) = [ ...
      -0.133333333333333333, ...
       0.075000000000000000, ...
     ];

  elseif ( rule == 4 )

    suborder_xyzz(1:4,1:suborder_num) = [ ...
      0.568430584196844400,  0.143856471934385200  ...
      0.143856471934385200,  0.143856471934385200; ...
      0.500000000000000000,  0.500000000000000000  ...
      0.000000000000000000,  0.000000000000000000  ...
    ]';

    suborder_w(1:suborder_num) = [ ...
      0.0362941783134009000,  ...
      0.00358165890217718333  ...
    ];

  elseif ( rule == 5 )

    suborder_xyzz(1:4,1:suborder_num) = [ ...
      0.250000000000000000,   0.250000000000000000   ...
      0.250000000000000000,   0.250000000000000000;  ...
      0.785714285714285714,   0.0714285714285714285  ...
      0.0714285714285714285,  0.0714285714285714285; ...
      0.399403576166799219,   0.399403576166799219   ...
      0.100596423833200785,   0.100596423833200785   ...
    ]';

    suborder_w(1:suborder_num) = [ ...
     -0.0131555555555555556,  ...
      0.00762222222222222222, ...
      0.0248888888888888889   ...
    ];

  elseif ( rule == 6 )

    suborder_xyzz(1:4,1:suborder_num) = [ ...
      0.500000000000000000,   0.500000000000000000  ...
      0.000000000000000000,   0.000000000000000000; ...
      0.698419704324386603,   0.100526765225204467  ...
      0.100526765225204467,   0.100526765225204467; ...
      0.0568813795204234229,  0.314372873493192195  ...
      0.314372873493192195,   0.314372873493192195  ...
    ]';

    suborder_w(1:suborder_num) = [ ...
      0.00317460317460317450,  ...
      0.0147649707904967828, ...
      0.0221397911142651221   ...
    ];

  elseif ( rule == 7 )

    suborder_xyzz(1:4,1:suborder_num) = [ ...
      0.250000000000000000,  0.250000000000000000   ...
      0.250000000000000000,  0.250000000000000000;  ...
      0.00000000000000000,   0.333333333333333333   ...
      0.333333333333333333,  0.333333333333333333;  ...
      0.727272727272727273,  0.0909090909090909091  ...
      0.0909090909090909091, 0.0909090909090909091; ...
      0.0665501535736642813, 0.0665501535736642813  ...
      0.433449846426335728,  0.433449846426335728   ...
    ]';

    suborder_w(1:suborder_num) = [ ...
      0.0302836780970891856,  ...
      0.00602678571428571597, ...
      0.0116452490860289742,  ...
      0.0109491415613864534   ...
    ];

  elseif ( rule == 8 )

    suborder_xyzz(1:4,1:suborder_num) = [ ...
      0.356191386222544953,  0.214602871259151684   ...
      0.214602871259151684,  0.214602871259151684;  ...
      0.877978124396165982,  0.0406739585346113397  ...
      0.0406739585346113397, 0.0406739585346113397; ...
      0.0329863295731730594, 0.322337890142275646   ...
      0.322337890142275646,  0.322337890142275646;  ...
      0.0636610018750175299, 0.0636610018750175299  ...
      0.269672331458315867,  0.603005664791649076   ...
    ]';

    suborder_w(1:suborder_num) = [ ...
      0.00665379170969464506, ...
      0.00167953517588677620, ...
      0.00922619692394239843, ...
      0.00803571428571428248  ...
    ];

  elseif ( rule == 9)

    suborder_xyzz(1:4,1:suborder_num) = [ ...
      0.250000000000000000,    0.250000000000000000   ...
      0.250000000000000000,    0.250000000000000000;  ...
      0.765360423009044044,    0.0782131923303186549  ...
      0.0782131923303186549,   0.0782131923303186549; ...
      0.634470350008286765,    0.121843216663904411   ...
      0.121843216663904411,    0.121843216663904411;  ...
      0.00238250666073834549,  0.332539164446420554   ...
      0.332539164446420554,    0.332539164446420554;  ...
      0.500000000000000000,    0.500000000000000000   ...
      0.00000000000000000,     0.00000000000000000;   ...
      0.100000000000000000,    0.100000000000000000   ...
      0.200000000000000000,    0.600000000000000000   ...
    ]';

    suborder_w(1:suborder_num) = [ ...
      0.0182642234661087939,   ...
      0.0105999415244141609,   ...
     -0.0625177401143299494,   ...
      0.00489142526307353653,  ...
      0.000970017636684296702, ...
      0.0275573192239850917    ...
    ];

  elseif ( rule == 10 )

    suborder_xyzz(1:4,1:suborder_num) = [ ...
      0.250000000000000000,  0.250000000000000000   ...
      0.250000000000000000,  0.250000000000000000;  ...
      0.617587190300082967,  0.127470936566639015   ...
      0.127470936566639015,  0.127470936566639015;  ...
      0.903763508822103123,  0.0320788303926322960  ...
      0.0320788303926322960, 0.0320788303926322960; ...
      0.0497770956432810185, 0.0497770956432810185  ...
      0.450222904356718978,  0.450222904356718978;  ...
      0.183730447398549945,  0.183730447398549945   ...
      0.316269552601450060,  0.316269552601450060;  ...
      0.231901089397150906,  0.231901089397150906   ...
      0.0229177878448171174, 0.513280033360881072;  ...
      0.0379700484718286102, 0.0379700484718286102  ...
      0.730313427807538396,  0.193746475248804382   ...
    ]';

    suborder_w(1:suborder_num) = [ ...
     -0.0393270066412926145,   ...
      0.00408131605934270525,  ...
      0.000658086773304341943, ...
      0.00438425882512284693,  ...
      0.0138300638425098166,   ...
      0.00424043742468372453,  ...
      0.00223873973961420164   ...
    ];

  else

    fprintf ( 1, '\n' );
    fprintf ( 1, 'KEAST_SUBRULE - Fatal error!\n' );
    fprintf ( 1, '  Illegal RULE = %d\n', rule );
    error ( 'KEAST_SUBRULE - Fatal error!' );

  end
%
%  Renormalize the weights so that they sum to 1.
%
  suborder_w(1:suborder_num) = 6.0 * suborder_w(1:suborder_num);

  return
end