#!/usr/bin/env texlua

-- $Id: complex.lua 26 2024-02-04 22:42:59Z reinhard $

-- Copyright (C) 2020 Reinhard Kotucha <reinhard.kotucha@gmx.de>
-- 
-- 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., 51 Franklin Street, Fifth Floor, Boston,
-- MA 02110-1301, USA.

module(..., package.seeall)

-- Textual representation of complex numbers, for example "3 + 5i".
function cfmt (c, digits)
  local fmt = string.format
  -- The optional parameter digits restricts the number of fractional
  -- digits to be displayed.
  local format = '%g'
  if digits then format = '%.'..digits..'g' end

  if type(c) == 'table' then -- complex number
    if Im(c) >= 0 then
      return fmt(format..' + '..format..'i', Re(c), Im(c))
    else
      return fmt(format..' - '..format..'i', Re(c), -Im(c))
    end
  elseif type(c) == 'number' then -- real number
    return fmt(format, c)
  else
    return c
  end
end

function cprint (c, digits) print(cfmt(c, digits)) end

function Re (c) return c[1] end

function Im (c) return c[2] end

function conj (c) return { Re(c), -Im(c) } end

function abs (c)
  if type(c) == 'number' then
    return math.abs(c)
  else
    return math.sqrt(Re(c)^2 + Im(c)^2) -- Pythagoras
  end
end

function add (c1, c2)
  if type(c1) == 'number' then c1 = { c1, 0 } end
  if type(c2) == 'number' then c2 = { c2, 0 } end
  return { Re(c1) + Re(c2), Im(c1) + Im(c2) }
end

function sub (c1, c2)
  if type(c1) == 'number' then c1 = { c1, 0 } end
  if type(c2) == 'number' then c2 = { c2, 0 } end
  return { Re(c1) - Re(c2), Im(c1) - Im(c2) }
end

function mul (c1, c2)
  if type(c1) == 'number' then c1 = { c1, 0 } end
  if type(c2) == 'number' then c2 = { c2, 0 } end
  return { Re(c1)*Re(c2) - Im(c1)*Im(c2), Re(c1)*Im(c2) + Im(c1)*Re(c2) }
end

function div (c1, c2)
  if type(c1) == 'number' then c1 = { c1, 0 } end
  if type(c2) == 'number' then c2 = { c2, 0 } end
  local denom = Re(c2)^2 + Im(c2)^2
  if denom <= 0 then error('complex.div: denominator <= 0') end 
  return { (Re(c1)*Re(c2) + Im(c1)*Im(c2))/denom,
           (Im(c1)*Re(c2) - Re(c1)*Im(c2))/denom } 
end


-- Local Variables:
--  mode: Lua
--  lua-indent-level: 2
--  indent-tabs-mode: nil
--  coding: utf-8-unix
-- End:
-- vim:set tabstop=2 expandtab: