Skip to content
This repository has been archived by the owner on Apr 24, 2023. It is now read-only.

Commit

Permalink
Update to julia 0.6
Browse files Browse the repository at this point in the history
Mostly broadcast stuff, and some changes forthe sake of Lint.jl
  • Loading branch information
mdavezac committed Jul 9, 2017
1 parent 634f6b1 commit e0559ee
Show file tree
Hide file tree
Showing 10 changed files with 150 additions and 142 deletions.
4 changes: 2 additions & 2 deletions REQUIRE
Original file line number Diff line number Diff line change
@@ -1,5 +1,5 @@
julia 0.5
Unitful 0.0.4
julia 0.6
Unitful
DataFrames
CoordinateTransformations
Lumberjack
Expand Down
47 changes: 23 additions & 24 deletions src/Gruber.jl
Original file line number Diff line number Diff line change
Expand Up @@ -6,12 +6,11 @@ using Crystals.Utilities: cell_parameters
using Crystals.Structures: Crystal, volume
using Crystals: Log
using Unitful: unit, ustrip, Quantity
const max_no_change = 2

function no_opt_change_test(new, last)
const m_new = 16e0 * new;
const diff = new - last;
const m_new_plus_diff = m_new + diff;
const difference = new - last;
const m_new_plus_diff = m_new + difference;
const m_new_plus_diff_minus_m_new = m_new_plus_diff - m_new;
m_new_plus_diff_minus_m_new != 0;
end
Expand All @@ -31,8 +30,8 @@ function def_test(args; tolerance=default_tolerance)
end

function def_gt_0(args...; tolerance=default_tolerance)
zero, positive = def_test(args; tolerance=tolerance)
positive == 3 || (zero == 0 && positive == 1)
z₀, positive = def_test(args; tolerance=tolerance)
positive == 3 || (z₀ == 0 && positive == 1)
end


Expand All @@ -53,7 +52,7 @@ function n3_action(params::Vector, rinv::Matrix; tolerance=default_tolerance)
const j = params[5] -tolerance ? -1 : 1
const k = params[6] -tolerance ? -1 : 1
rinv[:, :] = rinv * [i 0 0; 0 j 0; 0 0 k]
params[4:end] = abs(params[4:end])
params[4:end] = abs.(params[4:end])
end

function n4_action(params::Vector, rinv::Matrix; tolerance=default_tolerance)
Expand All @@ -73,31 +72,31 @@ function n4_action(params::Vector, rinv::Matrix; tolerance=default_tolerance)
end
end
rinv[:, :] = rinv * update
params[4:end] = -abs(params[4:end])
params[4:end] = -abs.(params[4:end])
end

function n5_action(params::Vector, rinv::Matrix; tolerance=default_tolerance)
const sign = params[4] > tolerance ? -1 : 1
rinv[:, :] = rinv * [1 0 0; 0 1 sign; 0 0 1]
params[3] += params[2] + sign * params[4]
params[4] += 2sign * params[2]
params[5] += sign * params[6]
const pos_or_neg = params[4] > tolerance ? -1 : 1
rinv[:, :] = rinv * [1 0 0; 0 1 pos_or_neg; 0 0 1]
params[3] += params[2] + pos_or_neg * params[4]
params[4] += 2pos_or_neg * params[2]
params[5] += pos_or_neg * params[6]
end

function n6_action(params::Vector, rinv::Matrix; tolerance=default_tolerance)
const sign = params[5] > tolerance ? -1 : 1
rinv[:, :] = rinv * [1 0 sign; 0 1 0; 0 0 1]
params[3] += params[1] + sign * params[5]
params[4] += sign * params[6]
params[5] += 2sign * params[1]
const pos_or_neg = params[5] > tolerance ? -1 : 1
rinv[:, :] = rinv * [1 0 pos_or_neg; 0 1 0; 0 0 1]
params[3] += params[1] + pos_or_neg * params[5]
params[4] += pos_or_neg * params[6]
params[5] += 2pos_or_neg * params[1]
end

function n7_action(params::Vector, rinv::Matrix; tolerance=default_tolerance)
const sign = params[6] > tolerance ? -1 : 1
rinv[:, :] = rinv * [1 sign 0; 0 1 0; 0 0 1]
params[2] += params[1] + sign * params[6]
params[4] += sign * params[5]
params[6] += 2sign * params[1]
const pos_or_neg = params[6] > tolerance ? -1 : 1
rinv[:, :] = rinv * [1 pos_or_neg 0; 0 1 0; 0 0 1]
params[2] += params[1] + pos_or_neg * params[6]
params[4] += pos_or_neg * params[5]
params[6] += 2pos_or_neg * params[1]
end

function n8_action(params::Vector, rinv::Matrix)
Expand Down Expand Up @@ -145,7 +144,7 @@ function gruber{T <: Number}(cell::Matrix{T};
const metric = transpose(cell) * cell
params = vcat(diag(metric), [2metric[2, 3], 2metric[1, 3], 2metric[1, 2]])
rinv = eye(size(metric, 1))
nochange, previous = 0, -params[1:3]
no_change, previous = 0, -params[1:3]
iteration::Int = 0
for iteration in 1:itermax
condition0 =
Expand All @@ -160,7 +159,7 @@ function gruber{T <: Number}(cell::Matrix{T};
n3_action(params, rinv; tolerance=ε)
else
n4_action(params, rinv; tolerance=ε)
if all(abs(previous - params[1:3]) .< ε)
if all(abs.(previous .- params[1:3]) .< ε)
no_change += 1
else
no_change = 0
Expand Down
28 changes: 14 additions & 14 deletions src/SNF.jl
Original file line number Diff line number Diff line change
Expand Up @@ -3,16 +3,16 @@ using Unitful: Quantity, ustrip
export smith_normal_form

function choose_pivot!{T <: Integer}(left::Matrix{T}, smith::Matrix{T},
step::Integer, index::Integer=1)
istep::Integer, index::Integer=1)
index > size(smith, 2) && return 0
while all(smith[:, index] .== 0)
index += 1
index size(smith, 2) || return 0;
end
if smith[step, index] == 0
if smith[istep, index] == 0
k = findfirst(x -> x 0, smith[:, index])
left[step, :], left[k, :] = deepcopy(left[k, :]), deepcopy(left[step, :])
smith[step, :], smith[k, :] = deepcopy(smith[k, :]), deepcopy(smith[step, :])
left[istep, :], left[k, :] = deepcopy(left[k, :]), deepcopy(left[istep, :])
smith[istep, :], smith[k, :] = deepcopy(smith[k, :]), deepcopy(smith[istep, :])
end
index
end
Expand All @@ -25,11 +25,11 @@ function improve_col_pivot!{T <: Integer}(left::Matrix{T}, smith::Matrix{T},

β, σ, τ = gcdx(smith[row, col], smith[k, col])
α = smith[row, col] / β
γ = smith[k, col] / β
γₒ = smith[k, col] / β

Lp = eye(T, size(left, 1))
Lp[row, [row, k]] = [σ, τ]
Lp[k, [row, k]] = [-γ, α]
Lp[k, [row, k]] = [-γₒ, α]

left[:] = Lp * left
smith[:] = Lp * smith
Expand All @@ -44,11 +44,11 @@ function improve_row_pivot!{T <: Integer}(smith::Matrix{T}, right::Matrix{T},

β, σ, τ = gcdx(smith[row, col], smith[row, k])
α = smith[row, col] / β
γ = smith[row, k] / β
γ = smith[row, k] / β

Rp = eye(T, size(right, 1))
Rp[[col, k], col] = [σ, τ]
Rp[[col, k], k] = [-γ, α]
Rp[[col, k], k] = [-γ, α]

right[:] = right * Rp
smith[:] = smith * Rp
Expand Down Expand Up @@ -81,14 +81,14 @@ end

function diagonalize_all_entries!{T <: Integer}(left::Matrix{T}, smith::Matrix{T}, right::Matrix{T})
@assert size(smith, 1) == size(smith, 2)
step = 1
col = choose_pivot!(left, smith, step)
istep = 1
col = choose_pivot!(left, smith, istep)
while 0 < col size(smith, 2)
diagonalize_at_entry!(left, smith, right, step, col)
diagonalize_at_entry!(left, smith, right, istep, col)

step += 1
istep += 1
col += 1
col = choose_pivot!(left, smith, step, col)
col = choose_pivot!(left, smith, istep, col)
end
return smith, left, right
end
Expand Down Expand Up @@ -126,7 +126,7 @@ function smith_normal_form{T <: Integer}(matrix::Matrix{T})
move_zero_entries!(smith, right)
make_divisible!(left, smith, right)
left = convert(Matrix{T}, diagm([u 0 ? 1: -1 for u in diag(smith)])) * left
smith = abs(smith)
smith = abs.(smith)
smith, left, right
end

Expand Down
7 changes: 3 additions & 4 deletions src/SpaceGroup.jl
Original file line number Diff line number Diff line change
Expand Up @@ -55,7 +55,6 @@ Implementation taken from [ENUM](http://enum.sourceforge.net/).
function point_group{T <: Number}(cell::AbstractMatrix{T};
tolerance::Real=default_tolerance)
@assert size(cell, 1) == size(cell, 2)
const ndims = size(cell, 1)

avecs, bvecs, cvecs = potential_equivalents(cell, tolerance=tolerance)
result = Matrix{T}[eye(cell)]
Expand Down Expand Up @@ -112,7 +111,7 @@ function inner_translations_impl(fractional::AbstractMatrix,
species[site] atom_type && continue

translation = into_voronoi(fractional[:, site] - atom_center, grubcell)
all(abs(translation) .< tolerance) && continue
all(abs.(translation) .< tolerance) && continue

is_mapping = true
for mapping eachindex(species)
Expand Down Expand Up @@ -263,7 +262,7 @@ function primitive_impl(cell::AbstractMatrix,
det(ustrip(trial)) < 0e0 && Log.error("Negative volume")

int_cell = inv(trial) * cell
all(abs(int_cell - round(Integer, int_cell) .< 1e-8)) || continue
all(abs.(int_cell - round.(Integer, int_cell) .< 1e-8)) || continue

new_cell = trial
V = volume(trial)
Expand Down Expand Up @@ -344,7 +343,7 @@ function space_group_impl(cell::AbstractMatrix,
Log.info(
"$(length(result)) symmetry operations found, with " *
"$(count(result) do op
all(abs(ustrip(op(zeros(eltype(cartesian), size(translations, 1))))) .< 1e-8)
all(abs.(ustrip(op(zeros(eltype(cartesian), size(translations, 1))))) .< 1e-8)
end) " * "pure symmetries."
)
[u for u in result]
Expand Down
41 changes: 22 additions & 19 deletions src/Structures.jl
Original file line number Diff line number Diff line change
Expand Up @@ -16,9 +16,9 @@ using Unitful: dimension, unit
const RESERVED_COLUMNS = [:position, :fractional, :cartesian, :x, :y, :z]

""" Top type node for Crystals """
abstract AbstractCrystal
abstract type AbstractCrystal end

typealias RowIndices{T <: Integer} Union{T, AbstractVector{T}, Range{T}, Colon}
RowIndices{T <: Integer} = Union{T, AbstractVector{T}, Range{T}, Colon}

""" Describe a crystalline structure """
type Crystal{T <: Number, D, U, P <: Number} <: AbstractCrystal
Expand Down Expand Up @@ -65,8 +65,8 @@ function Crystal{C, D, U, P <: Number}(cell::Matrix{Quantity{C, D, U}},
end

function Crystal{C, D, U}(cell::Matrix{Quantity{C, D, U}}; kwargs...)
const dpositions = filter(x -> x[1] (:position, :positions), kwargs)
const tpositions = filter(x -> x[1] (:tposition, :tpositions), kwargs)
const dpositions = collect(filter(x -> x[1] (:position, :positions), kwargs))
const tpositions = collect(filter(x -> x[1] (:tposition, :tpositions), kwargs))
if length(dpositions) 0 && length(tpositions) 0
const positions = hcat((x[2] for x in dpositions)...,
transpose(hcat([x[2] for x in tpositions]...)))
Expand All @@ -75,12 +75,11 @@ function Crystal{C, D, U}(cell::Matrix{Quantity{C, D, U}}; kwargs...)
elseif length(tpositions) 0
const positions = transpose(hcat([x[2] for x in tpositions]...))
else
const positions = reshape(Matrix{Quantity{C, D, U}}(),
(size(cell, 1), 0))
const positions = Matrix{Quantity{C, D, U}}(size(cell, 1), 0)
end

filter!(x -> x[1] (:position, :positions, :tposition, :tpositions), kwargs)
Crystal(cell, positions; kwargs...)
leftover = filter(x -> x[1] (:position, :positions, :tposition, :tpositions), kwargs)
Crystal(cell, positions; leftover...)
end

""" Underlying physical units of the crystal """
Expand Down Expand Up @@ -128,17 +127,17 @@ volume(crystal::Crystal) = volume(crystal.cell)

function are_compatible_lattices(cell_a::Matrix, cell_b::Matrix;
digits=12, rtol=default_tolerance, atol=default_tolerance)
isinteger(round(inv(cell_a) * cell_b, digits)) &&
all(isinteger.(round.(inv(cell_a) * cell_b, digits))) &&
isapprox(ustrip(volume(cell_a)), ustrip(volume(cell_b)); rtol=rtol, atol=atol)
end
function are_compatible_lattices(crystal::Crystal, cell::Matrix; kwargs...)
are_compatible_lattices(crystal.cell, cell)
are_compatible_lattices(crystal.cell, cell; kwargs...)
end
function are_compatible_lattices(crystal_a::Crystal, crystal_b::Crystal; kwargs...)
are_compatible_lattices(crystal_a.cell, crystal_b.cell)
are_compatible_lattices(crystal_a.cell, crystal_b.cell; kwargs...)
end
function are_compatible_lattices(cell::Matrix, crystal::Crystal; kwargs...)
are_compatible_lattices(cell, crystal.cell)
are_compatible_lattices(cell, crystal.cell; kwargs...)
end
"""
are_compatible_lattices(lattices...)
Expand Down Expand Up @@ -201,7 +200,7 @@ end
position_for_crystal(::Val{:fractional}, ::AbstractMatrix, pos::AbstractArray) = pos
position_for_crystal(::Val{:cartesian}, cell::AbstractMatrix, p::AbstractArray) = cell * p
function position_for_crystal{T <: Quantity}(::Val{:cartesian},
cell::AbstractMatrix,
::AbstractMatrix,
positions::AbstractArray{T})
positions
end
Expand Down Expand Up @@ -383,7 +382,7 @@ function Base.getindex(crystal::Crystal, ::Colon, symbol::Symbol)
end

Base.getindex(crystal::Crystal, ::Colon, ::Colon) = copy(crystal)
Base.getindex(crystal::Crystal, row::Any, col::Colon) = Base.getindex(crystal, row)
Base.getindex(crystal::Crystal, row::Any, ::Colon) = Base.getindex(crystal, row)

function Base.getindex(crystal::Crystal, index::RowIndices, symbols::AbstractVector{Symbol})
const specials = symbols RESERVED_COLUMNS
Expand Down Expand Up @@ -428,7 +427,7 @@ end
Base.names(crystal::Crystal) = push!(names(crystal.properties), :position)
Base.size(crystal::Crystal) = (r = size(crystal.properties); (r[1], r[2] + 1))
Base.size(crystal::Crystal, i::Integer) = size(crystal)[i]
Base.ndims(crystal::Crystal) = 2
Base.ndims(::Crystal) = 2
DataFrames.nrow(crystal::Crystal) = size(crystal.positions, 2)
DataFrames.ncol(crystal::Crystal) = ncol(crystal.properties) + 1
Base.endof(crystal::Crystal) = nrow(crystal)
Expand Down Expand Up @@ -588,6 +587,10 @@ function DataFrames.deleterows!(crystal::Crystal, row::Integer)
crystal.positions = crystal.positions[:, rows]
end

function DataFrames.deleterows!(crystal::Crystal, ::Colon)
deleterows!(crystal, 1:length(crystal))
end

function DataFrames.deleterows!(crystal::Crystal, rows::RowIndices)
deleterows!(crystal.properties, rows)
prows = filter(i -> i rows, 1:nrow(crystal))
Expand Down Expand Up @@ -644,15 +647,15 @@ Rounds the cell and positions of a crystal. See `round` for possible parameters.
"""
function round!{T, D, U, TT, DD, UU}(crystal::Crystal{T, D, U, Quantity{TT, DD, UU}},
args...)
crystal.cell = round(reinterpret(T, crystal.cell), args...) * unit(Quantity{T, D, U})
crystal.cell = round.(reinterpret(T, crystal.cell), args...) * unit(Quantity{T, D, U})
const punit = unit(Quantity{TT, DD, UU})
crystal[:position] = round(reinterpret(TT, crystal[:position]), args...) * punit
crystal[:position] = round.(reinterpret(TT, crystal[:position]), args...) * punit
crystal
end

function round!{T, D, U, TT}(crystal::Crystal{T, D, U, TT}, args...)
crystal.cell = round(reinterpret(T, crystal.cell), args...) * unit(Quantity{T, D, U})
crystal[:position] = round(crystal[:position], args...)
crystal.cell = round.(reinterpret(T, crystal.cell), args...) * unit(Quantity{T, D, U})
crystal[:position] = round.(crystal[:position], args...)
crystal
end

Expand Down
Loading

0 comments on commit e0559ee

Please sign in to comment.