conjugate gradient method out of memory

[Yoshi Rokuko]

can someone see why i run out of memory? I know slices might
be wrong here, but the method works the first 90 iterations
and all slices should be reused. I think I missunderstood
something here:

for i := 0; i < sweeps; i++ {
s = laplace.Operator(shift, mass, nn)
prodDot = DotProd(shift, s)
shiftScale = dotProd/prodDot
for k, sk := range s {
phi[k] += shiftScale * shift[k]
rho[k] -= shiftScale * sk
}
prodDot = DotProd(rho, rho)
fmt.Println(i, math.Sqrt(prodDot/rel))
if math.Sqrt(prodDot/rel) < tol {
return phi
}
scaleShift = prodDot/dotProd
for k, rhok := range rho {
shift[k] = rhok + scaleShift * shift[k]
}
dotProd = prodDot
}
(full code see at the end)

my first version was using s:=laplace.. and prodDot:=DotProd..
and so on, but that's about the same ...

Thank you, best regards, Yoshi

package main

import (
"fmt"
"./gridgeom"
"./laplace"
"math"
"os"
"rand"
)

func main() {
ndim := 4
lsize := make([]int, ndim + 1)
for i, _ := range lsize {
lsize[i] = 31
}
nn, cols, nvol := gridgeom.PeriodicBound(lsize, ndim)
if len(nn) < ndim {
fmt.Fprintf(os.Stderr, "error: gridgeom.PeriodicBound() no next neighbor field\n")
os.Exit(1)
}
// Field initialisation
rho := make([]float64, nvol)
for i := 0; i < nvol; i++ {
if nn[i*cols] == 1 {
rho[i] = rand.NormFloat64()
}
}
fmt.Printf("set terminal png\nset output 'convergence.png'\nset title '%v%s %v %s'\nset logscale y\nplot '-'\n", ndim, "D grid with", nvol, "points")
phi := ConjGrad(rho, nn, 0, 10E-14)
if len(phi) < ndim {
fmt.Fprintf(os.Stderr, "error: lineq.ConjGrad() did not converge\n")
os.Exit(1)
}
}

func ConjGrad(rho []float64, nn []int, mass, tol float64) []float64 {
var prodDot, shiftScale, scaleShift float64
sweeps := 10000
nvol := len(rho)
s := make([]float64, nvol)
shift := make([]float64, nvol)
phi := make([]float64, nvol)
copy(shift, rho)
dotProd := DotProd(rho, rho)
rel := dotProd
for i := 0; i < sweeps; i++ {
s = laplace.Operator(shift, mass, nn)
prodDot = DotProd(shift, s)
shiftScale = dotProd/prodDot
for k, sk := range s {
phi[k] += shiftScale * shift[k]
rho[k] -= shiftScale * sk
}
prodDot = DotProd(rho, rho)
fmt.Println(i, math.Sqrt(prodDot/rel))
if math.Sqrt(prodDot/rel) < tol {
return phi
}
scaleShift = prodDot/dotProd
for k, rhok := range rho {
shift[k] = rhok + scaleShift * shift[k]
}
dotProd = prodDot
}
return make([]float64, 1)
}

func DotProd(x, y []float64) (res float64) {
for i, xi := range x {
res += xi*y[i]
}
return
}

================================================
package gridgeom

func PeriodicBound(lsize []int, ndim int) ([]int, int, int) {
if len(lsize) != ndim + 1 {
return make([]int, 0), 0, 0
}
ibase := make([]int, ndim + 2)
ibase[1] = 1
for k := 1; k <= ndim; k++ {
ibase[k + 1] = ibase[k]*lsize[k]
}
nvol := ibase[ndim + 1]
ix := make([]int, ndim + 1)
for i, _ := range ix {
ix[i] = 0
}
cols := ndim*2 + 1
nn := make([]int, nvol*cols)
for i := 0; i < nvol; i++ {
nn[i*cols] = 1
}
for i := 0; i < nvol; i++ {
irow := i*cols
for k := 1; k <= ndim; k++ {
nn[irow + k] = i + ibase[k]
if ix[k] == (lsize[k] - 1) {
nn[irow + k] -= ibase[k + 1]
nn[irow] = 0
}
nn[irow + k + ndim] = i - ibase[k]
if ix[k] == 0 {
nn[irow + k + ndim] += ibase[k + 1]
nn[irow] = 0
}
}
for k := 1; k <= ndim; k++ {
ix[k]++
if ix[k] < lsize[k] {
break
}
ix[k] = 0
}
}
return nn, cols, nvol
}

======================================================
package laplace

func Operator(phi []float64, mass float64, nn []int) []float64 {
nvol := len(phi)
cols := len(nn)/nvol
phiout := make([]float64, nvol)
scale := mass * mass + float64(cols - 1)
for i := 0; i < nvol; i++ {
if nn[i*cols] == 0 {
// Dirichlet
phiout[i] = 0
} else {
phiout[i] = phi[i] * scale
for k := 1; k < cols; k++ {
phiout[i] -= phi[nn[i*cols + k]]
}
}
}
return phiout
}

[Olreich]

I'm not really sure, but I think it's your laplace.Operator that's
running you out of memory. this could be because the GC isn't cleaning
it up well, or just that it's creating 10000 [279]float64 which is way
more than 2GB allows

If it is Operator creating a new [279]float64 in phiout, here's a
proposed solution that would just rewrite phi and not have any new
float64 arrays running around:

main.47: s = laplace.Operator(&shift, mass, nn)

laplace:

func Operator(phi *[]float64, mass float64, nn []int) []float64 {

nvol := len(phi)
cols := len(nn)/nvol

temp := float64(0.0)

scale := mass * mass + float64(cols - 1)
for i := 0; i < nvol; i++ {
if nn[i*cols] == 0 {
// Dirichlet

temp = 0
} else {
temp = phi[i] * scale

for k := 1; k < cols; k++ {

temp -= phi[nn[i*cols + k]]
}
}
phi[i] = temp;
}
return phi
}

Hope that helps.

[peterGo]

Yoshi,

I ran our code ran successfully with 173 iterations, on a machine with
only 500MB of real memory running Ubuntu 9.10 32-bit.

170 1.6710767062277957e-13
171 1.3244048766755584e-13
172 1.0522709208641284e-13
173 8.406297772841032e-14

What is your GOOS and GOARCH? What does Mercurial hg id show for your
GOROOT? How much real memory do you have?

Peter

[peterGo]

> I ran our code ran successfully with 173 iterations, on a machine with
> only 500MB of real memory running Ubuntu 9.10 32-bit.

S/B: I ran your code successfully, with 173 iterations, on a machine

[peterGo]

Olreich

> it's creating 10000 [279]float64 which is way
> more than 2GB allows

float64 is 64 bits or 8 bytes (64 / 8). Therefore, 10, 000 * 279 *
(64 / 8) = 22,320,000 bytes or 22.3MB.

Peter

[peterGo]

Yoshi,

These are some results from runtime.MemStats.

386: Alloc 1,355,484,464; TotalAlloc 1,359,973,936
amd64: Alloc 93,819,168; TotalAlloc 1,351,979,240

Peter

[Yoshi Rokuko]

;::::::::::::::::::::::: peterGo :
> Yoshi,
>
> I ran your code successfully with 173 iterations, on a machine with

> only 500MB of real memory running Ubuntu 9.10 32-bit.
>
> 170 1.6710767062277957e-13
> 171 1.3244048766755584e-13
> 172 1.0522709208641284e-13
> 173 8.406297772841032e-14
>
> What is your GOOS and GOARCH? What does Mercurial hg id show for your
> GOROOT? How much real memory do you have?
>
> Peter
>

Strange, I have 500MB memory (something like 190MB free) on slackware 13.0 32bit
It stops at something like 90 1e-7.

I use:
export GOROOT=$HOME/go
export GOARCH=386
export GOOS=linux
export GOBIN=$HOME/755

hg id: f776656df34c

Best regards, Yoshi

[Yoshi Rokuko]

Thank you, I will test that.

I think my C version works a bit like that, I am still unsure about
that slice thing ...

Best regards, Yoshi

;::::::::::::::::::::::: Olreich :

[Olreich]

xD I skipped kilobytes in my calc, nvm! *hides*

[peterGo]

Yoshi,

I've constructed a simple test case, which seems to show that the 6g
code, for GOARCH=amd64, regularly runs the garbage collector, while 8g
code, forGOARCH=386, does not, which would be a bug.

package main

import (
"fmt"

"runtime"
)

func fs() []float64 {
r := make([]float64, 923521)
return r
}

func main() {
s := fs()
for i := 0; i < 100; i++ {
s = fs()
m := runtime.MemStats
fmt.Printf("i %d; Alloc %d; TotalAlloc %d\n", i, m.Alloc,
m.TotalAlloc)
}
_ = s
}


This would explain the final Alloc size disparity in my earlier
results:

386: Alloc 1,355,484,464; TotalAlloc 1,359,973,936
amd64: Alloc 93,819,168; TotalAlloc 1,351,979,240

I have some more checking to do before filing a bug report. I don't
see anything obvious in the assembler code.

Peter

[peterGo]

Yoshi,

I filed the results of my tests as: Issue 909: GOARCH=386 Garbage
Collection Is Ineffectual: mmap: errno=0xc
http://code.google.com/p/go/issues/detail?id=909

Peter

[Yoshi Rokuko]

> Yoshi,
>
> I filed the results of my tests as: Issue 909: GOARCH=386 Garbage
> Collection Is Ineffectual: mmap: errno=0xc
> http://code.google.com/p/go/issues/detail?id=909
>
> Peter

Thank you very much Peter!

Best regards, Yoshi

[andrey mirtchovski]

note that if you throw in a runtime.GC() call inside the for loop then
the code behaves similarly to the 6g one:

original code:

i 569; Alloc 4240432608; TotalAlloc 4408242120
i 570; Alloc 4247858704; TotalAlloc 4415960056
i 571; Alloc 4255284800; TotalAlloc 4423677992
mmap: errno=0xc

with runtime.GC() inside the for loop:

i 996; Alloc 289776696; TotalAlloc 7668478936
i 997; Alloc 289776696; TotalAlloc 7676160008
i 998; Alloc 289776696; TotalAlloc 7683841080
i 999; Alloc 289776696; TotalAlloc 7691522152
end; Alloc 282387512; TotalAlloc 7691814040

[peterGo]

Andrey,

Exactly, the garbage collector isn't kicking in properly.

I believe the memory manager uses several storage pools, ranked by
allocation size, with different allocation and management mechanisms.
If you reduce the allocation size from 923521 to 92352 or to 9235, the
garbage collector is effective again. Large memory chunks seem to be
the problem, which is probably why people haven't encountered the
problem before.

Peter