Problem with integration of fields

Hi,
I'm having a problem with the integration of the electric field. I'm
trying to calculate the effective area of an optical fiber. This is
defined as:

( integral ( abs(E)^2 ) )^2
----------------------------------
  integral ( abs(E)^4 )

where the integrals are over the area in the x-y plane and the E
fields are the z component.

When I do this I get values which are wrong by about a factor of
2-2.5. The following script shows this for a very simple example of a
silica strand (1.25 micrometers in radius) surrounded by air. For this
case the analytical solution can be found. For a wavelength of 1.06
micrometers, the propagation constant beta (z component of k) is 8.41
inverse micrometers and the effective area is 1.6 square micrometers.

With the MPB script below I get very good agreement for beta, but the
effective area is wrong. The output given is:

DATA 1 8.41113129877697 3.8272694189535

Any help would be very much appreciated.
Best regards,
John Travers

(define-param strand_radius 1.25) ; strand radius [micrometers]
(define-param wl 1.06)            ; wavelength [micrometers]
; setup computation
(set! geometry-lattice (make lattice (size 6 6 no-size)))
(set! resolution (vector3 32 32))
(set! tolerance 1e-8)
; structure
(define silica (make dielectric (epsilon 2.102)))
(set! default-material air)
(set! geometry (list (make cylinder (center 0 0 0) (radius strand_radius)
                                    (height infinity) (material silica))))
; variable definitions
(define pi (acos -1))
(define area 0)
; band function to get effective area
(define (proc_bands i)
  (get-efield i)
  (set! area ( * (/ (expt (compute-field-integral
                            (lambda (F eps r)
                              (expt (magnitude (vector3-z F)) 2))) 2)
                    (compute-field-integral
                      (lambda (F eps r)
                        (expt (magnitude (vector3-z F)) 4))))))
  (print "\nDATA " i " " (* 2.0 pi (vector3-z current-k)) " " area "\n"))
(find-k NO-PARITY (/ 1 wl) 1 1 (vector3 0 0 1) 1e-4 8.0 1.0 20.0 proc_bands)