#!/usr/bin/python -Wall

# ================================================================
# John Kerl
# kerl.john.r@gmail.com
# 2007-01-12
# ================================================================
# This software is released under the terms of the GNU GPL.
# Please see LICENSE.txt in the same directory as this file.
# ================================================================

# ================================================================
# John Kerl
# kerl.john.r@gmail.com
# 2007-01-12
#
# Note 2007-05-08:  The most recent AMS Notices has a nice article about better
# ways to do this -- namely, doing it uniformly distributed.  The matrices
# generated by this program are "random" in the sense of "different stuff each
# time", but it's not clear what the density function is.
# ================================================================

from __future__ import division # 1/2 = 0.5, not 0.
import sys
import random
import math
from sackmat_m import *

# ----------------------------------------------------------------
def randscalar():
	#lo    = -4.0
	#range =  8.0 # -10 to 10.
	#return lo + random.random() * range
	return random.gauss(0.0, 1.0)

# ----------------------------------------------------------------
def randmat(m, n):
	A = make_zero_matrix(m, n)
	for i in range(0, m):
		for j in range(0, n):
			A[i][j] = randscalar()
	return A

# ----------------------------------------------------------------
def randsqmat(n):
	return randmat(n, n)

# ----------------------------------------------------------------
def randgl(n):
	tol = 1e-10
	k = 1
	while (1):
		A = randmat(n, n)
		d = A.det()
		if (abs(d) >= tol):
			break
	return A

# ----------------------------------------------------------------
def randsl(n):
	A = randgl(n)
	d = A.det()
	absd = abs(d)
	absdroot = absd ** (1.0/n)
	for i in range(0, n):
		for j in range(0, n):
			A[i][j] /= absdroot
	if (d < 0.0):
		for j in range(0, n):
			A[0][j] = -A[0][j] 
	return A

# ----------------------------------------------------------------
def rando(n):
	A = randmat(n, n)
	A = gram_schmidt(A)
	return A

# ----------------------------------------------------------------
def randso(n):
	A = rando(n)
	d = A.det()
	if (d < 0.0):
		for j in range(0, n):
			A[0][j] = -A[0][j] 
	return A

# ----------------------------------------------------------------
def randnonneg(m,n):
	A = make_zero_matrix(m, n)
	for i in range(0, m):
		for j in range(0, n):
			A[i][j] = random.random()
	return A

# ----------------------------------------------------------------
def randsymm(n):
	A = randmat(n, n)
	A = A + A.transpose()
	return A

# ----------------------------------------------------------------
def randsksymm(n):
	A = randmat(n, n)
	A = A - A.transpose()
	return A