File: gsl_assist.ml

package info (click to toggle)
orpie 1.5.2-2
  • links: PTS
  • area: main
  • in suites: bullseye, buster, sid
  • size: 1,924 kB
  • ctags: 2,720
  • sloc: ml: 13,872; ansic: 3,754; makefile: 310; sh: 11; python: 11
file content (70 lines) | stat: -rw-r--r-- 2,289 bytes parent folder | download | duplicates (2)
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
(*  Orpie -- a fullscreen RPN calculator for the console
 *  Copyright (C) 2003-2004, 2005, 2006-2007, 2010 Paul Pelzl
 *
 *  This program is free software; you can redistribute it and/or modify
 *  it under the terms of the GNU General Public License, Version 2,
 *  as published by the Free Software Foundation.
 *
 *  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
 *
 *  Please send bug reports, patches, etc. to Paul Pelzl at 
 *  <pelzlpj@gmail.com>.
 *)

(* gsl_assist.ml
 *
 * Contains a bunch of helper functions that make ocamlgsl a bit easier
 * to work with. *)


(* simplifying functions for use with ocamlgsl bindings *)
let cmpx_of_int i   = {Complex.re=Big_int.float_of_big_int i; Complex.im=0.0}
let cmpx_of_float f = {Complex.re=f; Complex.im=0.0}
let cmat_of_fmat fm =
   let rows, cols = Gsl_matrix.dims fm and
   f_array = Gsl_matrix.to_array fm in
   let c_array = Array.map cmpx_of_float f_array in
   Gsl_matrix_complex.of_array c_array rows cols



(* 1-norm of matrix *)
let one_norm mat =
   let n, m = Gsl_matrix.dims mat in
   let maxval = ref (-1.0) in
   let sum = ref 0.0 in
   for j = 0 to pred m do
      sum := 0.0;
      for i = 0 to pred n do
         sum := !sum +. (abs_float mat.{i, j})
      done;
      if !sum > !maxval then
         maxval := !sum
      else
         ()
   done;
   !maxval



(* solve a complex linear system using LU decomposition *)
let solve_complex_LU ?(protect=true) mat b =
  let mA = Gsl_vectmat.cmat_convert ~protect mat in
  let vB = (`CV (Gsl_vector_complex.copy b)) in
  let (len, _) = Gsl_vectmat.dims mA in
  let p = Gsl_permut.create len in
  let _ = Gsl_linalg.complex_LU_decomp mA p in
  let x = Gsl_vector_complex_flat.create len in
  Gsl_linalg.complex_LU_solve mA p vB (`CVF x);
  Gsl_vector_complex_flat.to_complex_array x



(* arch-tag: DO_NOT_CHANGE_a19e0df2-6d6b-4925-87eb-be2a2926ffbb *)