Module: CodeRunner::Gs2::GSLMatrices

Included in:
CodeRunner::Gs2
Defined in:
lib/gs2crmod/gsl_data.rb

Instance Method Summary collapse

Instance Method Details

#es_heat_flux_over_ky_over_kx_gsl_matrix(options) ⇒ Object 



856
857
858
859
860
861
862
863
864
865
866
867
868
869
870
871
872
873
874
875
876
877
 
# File 'lib/gs2crmod/gsl_data.rb', line 856

def es_heat_flux_over_ky_over_kx_gsl_matrix(options)
Dir.chdir(@directory) do
		raise "Heat flux spectrum makes no sense for single modes" if @grid_option == "single"
		options.convert_to_index(:t) if options[:t] or options[:t_element]
		options[:t_index] ||= list(:t).keys.max
		#es_heat_by_k index order (in Fortran) is kx, ky, t
		es_heat_narray = netcdf_file.var("es_heat_by_k").get('start' => [0, 0, 0, options[:t_index] - 1], 'end' => [-1, -1, 0, options[:t_index] - 1])
		es_heat_narray.reshape!(*es_heat_narray.shape.slice(0..1))
		
		gm =  es_heat_narray.to_gm.move_cols_from_box_order
		if options[:limit]
			for i in 0...gm.shape[0]
				for j in 0...gm.shape[1]
# 						j+= extra if 
					gm[i, j] = [[gm[i,j], (options[:limit][0] or gm[i,j])].max, (options[:limit][1] or gm[i,j])].min
# 						mat[i, j+extra] = gm[i,-j] unless j==0
				end
			end
		end
		return gm
end
end

#growth_rate_over_ky_over_kx_gsl_matrix(options) ⇒ Object 



848
849
850
851
 
# File 'lib/gs2crmod/gsl_data.rb', line 848

def growth_rate_over_ky_over_kx_gsl_matrix(options)
		array = @growth_rate_at_ky_at_kx.values.map{|h| h.values}
		return GSL::Matrix.alloc(array.flatten, array.size, array[0].size)
end

#phi0_over_x_over_y_gsl_matrix(options) ⇒ Object 



959
960
961
962
963
964
965
966
967
968
969
970
971
972
973
974
975
976
977
978
979
980
981
982
983
984
985
986
987
988
989
990
991
992
993
994
995
996
997
998
999
1000
1001
1002
1003
1004
1005
1006
1007
1008
1009
1010
1011
1012
1013
1014
1015
1016
1017
1018
1019
1020
1021
1022
1023
1024
1025
1026
1027
1028
1029
1030
1031
 
# File 'lib/gs2crmod/gsl_data.rb', line 959

def phi0_over_x_over_y_gsl_matrix(options)
Dir.chdir(@directory) do

		options.convert_to_index(:t) if options[:t] or options[:t_element]
		options[:t_index] ||= list(:t).keys.max
		#phi2_by_mode index order (in Fortran) is kx, ky, t
		phi_re_narray = netcdf_file.var("phi0").get('start' => [0, 0, 0, options[:t_index] - 1], 'end' => [0, -1, -1, options[:t_index] - 1])
# 			ep phi_re_narray.shape
		phi_re_narray.reshape!(*phi_re_narray.shape.slice(1..2))
		# The narray has index order ky, kx, but we want kx, ky for historical reasons, hence the transpose. 
		gm_re = phi_re_narray.to_gm
		phi_im_narray = netcdf_file.var("phi0").get('start' => [1, 0, 0, options[:t_index] - 1], 'end' => [1, -1, -1, options[:t_index] - 1])
		phi_im_narray.reshape!(*phi_im_narray.shape.slice(1..2))
		# The narray has index order ky, kx, but we want kx, ky for historical imasons, hence the transpose. 
		gm_im = phi_im_narray.to_gm
		gm = GSL::Matrix::Complex.re_im(gm_re, gm_im)
# 			ep gm.shape

		if options[:no_zonal]
			
			for i in 0...gm.shape[1]
				gm[0,i] = GSL::Complex.alloc([0,0])
			end
		end
		if xres = (options[:xres] or options[:x_resolution])
			mat = GSL::Matrix::Complex.calloc(gm.shape[0], xres)
			extra = ((xres - gm.shape[1])).floor
			for i in 0...gm.shape[0]
				for j in 0...((gm.shape[1] + 1) / 2 )
# 						j+= extra if 
					mat[i, j] = gm[i,j]
					mat[i, j+extra] = gm[i,-j] unless j==0
				end
			end
			gm = mat
			
			
# 				gm = mat.vertcat(gm).vertcat(mat)
		end
		if yres = (options[:yres] or options[:y_resolution])
			mat = GSL::Matrix::Complex.calloc(yres, gm.shape[1])
			extra = ((yres - gm.shape[0])).floor
			for i in 0...gm.shape[0]
				for j in 0...gm.shape[1]
# 						j+= extra if 
					mat[i, j] = gm[i,j]
# 						mat[i, j+extra] = gm[i,-j] unless j==0
				end
			end
			gm = mat
			
			
# 				gm = mat.vertcat(gm).vertcat(mat)
		end
# 			ep gm_re, gm_im
# 			re = GSL::Complex.alloc([1.0, 0.0])
# 			gm = GSL::Matrix::Complex.calloc(*gm_re.shape)
# 			gm = gm_re * re  + gm_im * GSL::Complex.alloc([0.0, 1.0])
# 			gm_re, gm_im = fourier_transform_gm_matrix_complex_rows(gm_re, gm_im)
		
		gm = gm.backward_cols_c2c(true).backward_rows_cc2r(true)  
		if options[:limit]
			for i in 0...gm.shape[0]
				for j in 0...gm.shape[1]
# 						j+= extra if 
					gm[i, j] = [[gm[i,j], options[:limit][0]].max, options[:limit][1]].min
# 						mat[i, j+extra] = gm[i,-j] unless j==0
				end
			end
		end
		return gm
end
end

#spectrum_over_ky_over_kpar_gsl_matrix(options) ⇒ Object 



912
913
914
915
916
917
918
919
920
921
922
923
924
925
926
927
928
929
930
931
932
933
934
935
936
937
938
939
940
941
942
943
944
945
946
947
948
949
950
951
952
953
954
955
956
957
 
# File 'lib/gs2crmod/gsl_data.rb', line 912

def spectrum_over_ky_over_kpar_gsl_matrix(options)
Dir.chdir(@directory) do

		#:re, :theta, :kx, :ky
		lkx = list(:kx)
	
		if options[:t_index] or options[:t]
			#extra option required is t_index
			raise CRFatal.new("write_phi_over_time is not enabled so this function won't work") unless @write_phi_over_time
			options.convert_to_index(self, :t)
		end
		temp = phi_av = (lkx.keys.map do |kx_index|		
			if options[:t_index]
				phi =  netcdf_file.var('phi_t').get({'start' => [0,0,kx_index-1,0, options[:t_index] - 1], 'end' => [-1,-1,kx_index-1,-1, options[:t_index] - 1]})
			else
				phi = netcdf_file.var('phi').get({'start' => [0, 0, kx_index - 1, 0], 'end' => [-1, -1, kx_index-1, -1]})
			end
			#ep phi.shape
			phi.reshape(*phi.shape.values_at(0,1,3))
		end).sum / lkx.size

		phi_t = phi_av.to_a #.map{|arr| arr.transpose}.transpose.map{|a| a.transpose}
		#ep 'phi_t', phi_t.size, phi_t[0].size, phi_t[0][0].size
		gvky = gsl_vector('ky')
		gm = GSL::Matrix.alloc(gvky.size, gsl_vector('theta').size)
		for ky_element in 0...gm.shape[0]
			#p phi_t[ky_element].transpose[0]
			spectrum = GSL::Vector::Complex.alloc(phi_t[ky_element]).forward.square
			if options[:no_kpar0]
				spectrum[0]=0.0
			end
			spectrum = spectrum.from_box_order
			#ep spectrum.shape
			spectrum = spectrum*gvky[ky_element]**2 unless options[:phi2_only]
			gm.set_row(ky_element, spectrum)
		end
		if options[:no_zonal]
			gm.row(0).set_all(0.0)
		end
		if options[:log]
			gm = gm.log
		end
		
		return gm
end
end

#spectrum_over_ky_over_kx_gsl_matrix(options) ⇒ Object 



878
879
880
881
882
883
884
885
886
887
888
889
890
891
892
893
894
895
896
897
898
899
900
901
902
903
904
905
906
907
908
909
910
911
 
# File 'lib/gs2crmod/gsl_data.rb', line 878

def spectrum_over_ky_over_kx_gsl_matrix(options)
Dir.chdir(@directory) do
		raise "Spectrum makes no sense for single modes" if @grid_option == "single"
		options.convert_to_index(:t) if options[:t] or options[:t_element]
		options[:t_index] ||= list(:t).keys.max
		#phi2_by_mode index order (in Fortran) is kx, ky, t
		phi_narray = netcdf_file.var("phi2_by_mode").get('start' => [0, 0, options[:t_index] - 1], 'end' => [-1, -1, options[:t_index] - 1])
		phi_narray.reshape!(*phi_narray.shape.slice(0..1))
		
		gm =  phi_narray.to_gm.move_cols_from_box_order
		if options[:times_kx4] or options[:times_kx2]
# 				puts 'normalising'
			vals = list(:kx).values.sort
			for i in 0...gm.shape[0]
				for j in 0...gm.shape[1]
# 						p vals[j]
					gm[i,j] =  gm[i,j] * (vals[j])**4 if options[:times_kx4]
					gm[i,j] =  gm[i,j] * (vals[j])**2 if options[:times_kx2]
				end
			end
		end
					if options[:no_zonal]
			
			for i in 0...gm.shape[1]
				gm[0,i] = 0.0
			end
		end
		if options[:log]
			gm = gm.log
		end

		return gm
end
end

#transient_amplification_over_ky_over_kx_gsl_matrix(options) ⇒ Object 



852
853
854
855
 
# File 'lib/gs2crmod/gsl_data.rb', line 852

def transient_amplification_over_ky_over_kx_gsl_matrix(options)
		array = @transient_amplification_at_ky_at_kx.values.map{|h| h.values}
		return GSL::Matrix.alloc(array.flatten, array.size, array[0].size)
end