# include <math.h>
# include <stdlib.h>

# include "traveling_wave_exact.h"

/******************************************************************************/

void traveling_wave_exact ( double x, double t, double *u, double *ut, 
  double *utt, double *ux, double *uxx )

/******************************************************************************/
/*
  Purpose:

    traveling_wave_exact() evaluates an exact solution of the 1D wave equation.

  Discussion:

    d^2u/dt^2 = c^2 d^2u/dx^2
    u[x,t] = a * sin ( k * x + w * t + phi )

  Licensing:

    This code is distributed under the MIT license. 

  Modified:

    11 April 2025

  Author:

    John Burkardt

  Input:

    real x, t: the position and time.

  Output:

    real u, ut, utt, ux, uxx: the values of the exact solution, its time
    derivative, and its first and second spatial derivatives at (x,t).
*/
{
  double a;
  double k;
  double phi;
  double t0;
  double tstop;
  double w;
  double xmax;
  double xmin;

  traveling_wave_parameters ( NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL,
    &a, &k, &w, &phi, &xmin, &xmax, &t0, &tstop );

  *u = a * sin ( k * x + w * t + phi );

  *ut = a * w * cos ( k * x + w * t + phi );

  *utt = - a * w * w * sin ( k * x + w * t + phi );

  *ux = a * k * cos ( k * x + w * t + phi );

  *uxx = - a * k * k * sin ( k * x + w * t + phi );

  return;
}
/******************************************************************************/

void traveling_wave_parameters ( double *a_in, double *k_in, double *w_in, 
  double *phi_in, double *xmin_in, double *xmax_in, double *t0_in, 
  double *tstop_in, double *a_out, double *k_out, double *w_out, 
  double *phi_out, double *xmin_out, double *xmax_out, double *t0_out, 
  double *tstop_out )

/******************************************************************************/
/*
  Purpose:

    traveling_wave_parameters() returns parameters for the traveling wave equation.

  Discussion:

    d^2u/dt^2 = c^2 d^2u/dx^2

    If input values are specified, this resets the default parameters.
    Otherwise, the output will be the current defaults.

  Licensing:

    This code is distributed under the MIT license.

  Modified:

    09 April 2025

  Author:

    John Burkardt

  Input:

    real a_in: the amplitude.

    real k_in: the wave number.

    real w_in: angular frequency.

    real phi_in: the phase angle.

    real xmin_in, xmax_in: the left and right interval endpoints.

    real t0_in: the initial time.

    real tstop_in: the final time.

  Output:

    real a_out: the amplitude.

    real k_out: the wave number.

    real w_out: angular frequency.

    real phi_out: the phase angle.

    real xmin_out, xmax_out: the left and right interval endpoints.

    real t0_out: the initial time.

    real tstop_out: the final time.
*/
{
  static double a_default = 2.0;
  static double k_default = 6.0;
  static double phi_default = 3.141592653589793/ 4.0;
  static double t0_default = 0.0;
  static double tstop_default = 5.0;
  static double w_default = 1.0;
  static double xmax_default = +1.0;
  static double xmin_default = -1.0;

/*
  New values, if supplied on input, overwrite the current values.
*/
  if ( a_in )
  {
    a_default = *a_in;
  }

  if ( k_in )
  {
    k_default = *k_in;
  }

  if ( w_in )
  {
    w_default = *w_in;
  }

  if ( phi_in )
  {
    phi_default = *phi_in;
  }

  if ( xmin_in )
  {
    xmin_default = *xmin_in;
  }

  if ( xmax_in )
  {
    xmax_default = *xmax_in;
  }

  if ( t0_in )
  {
    t0_default = *t0_in;
  }

  if ( tstop_in )
  {
    tstop_default = *tstop_in;
  }
/*
  Return values.
*/
  if ( a_out )
  {
    *a_out = a_default;
  }
  if ( k_out )
  {
    *k_out = k_default;
  }
  if ( w_out )
  {
    *w_out = w_default;
  }
  if ( phi_out )
  {
    *phi_out = phi_default;
  }
  if ( xmin_out )
  {
    *xmin_out = xmin_default;
  }
  if ( xmax_out )
  {
    *xmax_out = xmax_default;
  }
  if ( t0_out )
  {
    *t0_out = t0_default;
  }
  if ( tstop_out )
  {
    *tstop_out = tstop_default;
  }
  return;
}
/******************************************************************************/

void traveling_wave_residual ( double t, double x, double *r )

/******************************************************************************/
/*
  Purpose:

    traveling_wave_residual() computes the residual of the traveling wave equation.

  Discussion:

    d^2u/dt^2 = c^2 d^2u/dx^2

  Licensing:

    This code is distributed under the MIT license.

  Modified:

    11 April 2025

  Author:

    John Burkardt

  Input:

    real T, X: the time and position where the solution is evaluated.

  Output:

    double *r: the residual at that time and position.
*/
{
  double a;
  double c;
  double k;
  double phi;
  double t0;
  double tstop;
  double u;
  double ut;
  double utt;
  double ux;
  double uxx;
  double w;
  double xmax;
  double xmin;

  traveling_wave_parameters ( NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL,
    &a, &k, &w, &phi, &xmin, &xmax, &t0, &tstop );

  traveling_wave_exact ( t, x, &u, &ut, &utt, &ux, &uxx );

  c = w / k;

  *r = utt - c * c * uxx;

  return;
}